Deciphering the nonlocal entanglement entropy of fracton topological orders

TL;DR

This paper investigates the nonlocal entanglement entropy in fracton topological orders, deriving bounds that reveal geometry-dependent properties and potential distinctions from conventional topological phases.

Contribution

It introduces a lower bound for nonlocal entanglement entropy applicable to fracton models, extending understanding beyond traditional topological quantum field theory descriptions.

Findings

Lower bound applies to Abelian stabilizer models including fractons and topological orders.

Nonlocal entanglement entropy can be geometry-dependent and extensive in fracton models.

The lower bound remains stable under local perturbations.

Abstract

The ground states of topological orders condense extended objects and support topological excitations. This nontrivial property leads to nonzero topological entanglement entropy $S_{topo}$ for conventional topological orders. Fracton topological order is an exotic class of models which is beyond the description of TQFT. With some assumptions about the condensates and the topological excitations, we derive a lower bound of the nonlocal entanglement entropy $S_{nonlocal}$ (a generalization of $S_{topo}$). The lower bound applies to Abelian stabilizer models including conventional topological orders as well as type \Rom{1} and type \Rom{2} fracton models, and it could be used to distinguish them. For fracton models, the lower bound shows that $S_{nonlocal}$ could obtain geometry-dependent values, and $S_{nonlocal}$ is extensive for certain choices of subsystems, including some choices which always give zero for TQFT. The stability of the lower bound under local perturbations is discussed.

Figures

Equations

$$I(A:C|B)\equiv S_{AB}+S_{BC}-S_B-S_{ABC}$$

$$I(A:C|B){|}_{\rho }\ge S_{ABC}{|}_{\sigma }-S_{ABC}{|}_{\rho }$$

$$I(A:C|B){|}_{\rho }=(S_{BC}+S_{CD}-S_B-S_D){|}_{\rho }.$$

$$S_{nonlocal}\equiv (S_{BC}+S_{CD}-S_B-S_D){|}_{\rho }.$$

$$S_{nonlocal}\ge S_{ABC}{|}_{\sigma }-S_{ABC}{|}_{\rho }.$$

$$S_{nonlocal}=S_{ABC}{|}_{{\sigma }^{\ast }}-S_{ABC}{|}_{\rho }.$$

$$U_i|\psi ⟩=U_i^{def}|\psi ⟩.$$

$$U_i=U_{i\Omega }\otimes U_{i\bar{\Omega }}.$$

$$\bigg{(}\prod_{i=1}^{M}U_{i}^{k_{i}}\bigg{)}|\psi\rangle=e^{i\delta(\{k_{i}\})}\bigg{(}\prod_{i=1}^{M}[U_{i}^{def}]^{k_{i}}\bigg{)}|\psi\rangle,$$

$$W_i|\psi ⟩=|\psi ⟩i=1,\dots ,M.$$

$$W_i{U}_j=U_j{W}_i{e}^{i{\theta }_{ij}}\text{with}{\theta }_{ij}=\frac{2\pi }{n_i}{\delta }_{ij},$$

$$|\{k_i\};\psi ⟩\equiv \prod _{i=1}^M{U}_i^{k_i}|\psi ⟩\text{with}0\le k_i\le n_i-1$$

$$S_{ABC}{|}_{\sigma }-S_{ABC}{|}_{\rho }=-\sum _{I=1}^N{p}_I\ln p_I\le \ln N$$

$$S_{nonlocal}{|}_{\rho }\ge \ln N=\sum _{i=1}^M\ln n_i.$$

$$H=-\sum _s{A}_s-\sum _p{B}_p,$$

$$A_s=\prod _{r\in s}{X}_r{B}_p=\prod _{r\in p}{Z}_r,$$

$$U_1^2=U_2^2=W_1^2=W_2^2=1$$

$$W_i{U}_j=U_j{W}_i{e}^{i\pi {\delta }_{ij}}i,j=1,2.$$

$${\sigma }^{\ast }=\frac{1}{4}(\rho +U_1\rho U_1^{†}+U_2\rho U_2^{†}+U_1{U}_2\rho U_2^{†}{U}_1^{†})$$

$$H=-\sum _s{A}_s-\sum _p{B}_p;$$

$$A_s=\prod _{r\in s}{X}_r;B_p=\prod _{r\in p}{Z}_r.$$

$$S_{\Omega }=S_{\Omega ,local}+S_{\Omega ,topological}\Rightarrow S_{nonlocal}^{(\gamma )}=0.$$

$$S_{nonlocal}^{(\alpha )}=S_{nonlocal}^{(\beta )}=\ln 2;S_{nonlocal}^{(\gamma )}=0.$$

$$H=-\sum _c{A}_c-\sum _s(B_s^{(xy)}+B_s^{(yz)}+B_s^{(zx)}),$$

$$H=-\sum _{i,j,k}{h}_{(i,j,k)}^Z-\sum _{i,j,k}{h}_{(i,j,k)}^X,$$

$$h_{(0,0,0)}^Z$$

$$h_{(0,0,0)}^X$$

$$h^{Z}_{(0,0,0)}=Z\left(\begin{array}[]{c}1+y+xy\\ 1+z\end{array}\right);\quad h^{X}_{(0,0,0)}=X\left(\begin{array}[]{c}1+\bar{z}\\ 1+\bar{y}+\bar{x}\bar{y}\end{array}\right).$$

$$h^{Z}_{(0,0,0)}=Z\left(\begin{array}[]{c}1+f(x)y\\ 1+g(x)z\end{array}\right);\quad h^{X}_{(0,0,0)}=X\left(\begin{array}[]{c}1+\bar{g}(x)\bar{z}\\ 1+\bar{f}(x)\bar{y}\end{array}\right).$$

$$\big{(}S_{AB}-S_{A}\big{)}|_{\sigma^{\prime}}\simeq\big{(}S_{AB}-S_{A}\big{)}|_{\sigma},$$

$$\mathbf{U}\text{-}1,\mathbf{U}\text{-}3$$

Full text

I Introduction

Topological order Wen and Niu (1990) is a gapped quantum phase of matter beyond the description of the Landau-Ginzburg theory of symmetry breaking. Many of the early examples of topological orders (which we will refer to as conventional topological orders) share the following properties: robust ground state degeneracy which depends on the topology of the manifold Wen and Niu (1990), the ground states are locally indistinguishable Wen and Niu (1990); Bravyi et al. (2010); Bravyi and Hastings (2011), the existence of integer dimensional condensates and logical operators that can be topologically deformed Levin and Wen (2005), nontrivial braiding statistics of anyons or other topological excitations (or topologically charged excitations) e.g., excitations which could not be created alone by local operators Bravyi et al. (2011); Haah (2013), they are effectively described by topological quantum field theory (TQFT) at low temperatures Wen and Niu (1990), and they can be used to do fault-tolerant quantum information processing Kitaev (2003). And it is well known that, in 2D, a suitable linear combination of entanglement entropy with local contributions canceled is a topological invariance called the topological entanglement entropy Kitaev and Preskill (2006); Levin and Wen (2006). Topological entanglement entropy is a property of the ground state wave function and it has been used to identify quantum spin liquid phases Isakov et al. (2011). It also contains information about the ground state degeneracy Kim (2013) and the forms of low-energy excitations Kim and Brown (2015). Generalizations of topological entanglement entropy into 3D bulk Grover et al. (2011) and boundary Kim and Brown (2015) are studied.

On the other hand, there are recently discussed 3D exotic topological ordered models Chamon (2005); Bravyi et al. (2011); Haah (2011); Yoshida (2013); Vijay et al. (2015, 2016); Williamson (2016); Pretko (2017); Vijay (2017); Ma et al. (2017a); Hsieh and Halász (2017); Slagle and Kim (2017) that do not fit very well into the pictures above. These models have recently been classified into fracton topological orders Vijay et al. (2016). While fracton models have locally indistinguishable ground states when placed on nontrivial manifolds and the ground state degeneracy is robust under local perturbations Bravyi et al. (2010); Bravyi and Hastings (2011), the ground state degeneracy depends on the system size (geometry) rather than merely the topology of the manifold Bravyi et al. (2011); Haah (2013); Yoshida (2013). While fracton models possess topological excitations Bravyi et al. (2011); Haah (2013), these topological excitations are constrained to move in lower dimensional submanifolds rather than the whole system Vijay et al. (2015, 2016). The condensates and logical operators can be fractal dimensional Yoshida (2013); Vijay et al. (2016) instead of integer dimensional. These models are beyond the description of TQFT.

There are type I and type II fracton topological orders. The type I fracton models include the Chamon-Bravyi-Leemhuis-Terhal (CBLT) model Chamon (2005); Bravyi et al. (2011), the Majorana cubic model Vijay et al. (2015), and the X-cube model Vijay et al. (2016), etc.; they have integer dimensional condensates and logical operators. The type II fracton models include Haah’s code Haah (2011) and many of the fractal spin liquid models Yoshida (2013) (see Sec.III.4). Type II fracton models possess fractal condensates and logical operators and the excitations are fully immobile Vijay et al. (2016).

For the ground state entanglement properties of fracton models, a relation to the ground state degeneracy is implied in Kim (2013) and the entanglement renormalization group transformation of Haah’s code is studied in Haah (2014). In this work, we construct a direct analogy of topological entanglement entropy by doing linear combinations of entanglement entropies of different subsystems in such a way that the local contributions (from each boundary or corner of the subsystems) are canceled, and we call the linear combination $S_{nonlocal}$, the nonlocal entanglement entropy. While $S_{nonlocal}$ is topologically invariant for conventional topological orders, it is geometry-dependent for fracton models.

Explicitly, we choose a conditional mutual information form used by Kim and Brown Kim and Brown (2015) and define the nonlocal entanglement entropy $S_{nonlocal}\equiv(S_{BC}+S_{CD}-S_{C}-S_{D})|_{\rho}=I(A:C|B)|_{\rho}$. Where $\rho=|\psi\rangle\langle\psi|$ with $|\psi\rangle$ being a ground state. The whole system is the union of subsystems $A$, $B$, $C$, and $D$ with $A$ and $C$ separated by distance $l\gg\xi$, the correlation length. One can check that the local contributions from the boundaries are canceled, and this is why we use the name “nonlocal entanglement entropy.” This construction can be used in any dimensions and an example in 2D is shown in Fig.1. With several assumptions about the condensates and topological excitations, a lower bound of $S_{nonlocal}$ is derived.

When applied to known conventional Abelian topological orders e.g. the 2D toric code model and the 3D toric code model Dennis et al. (2002), the lower bound is topologically invariant and it is identical to the exact result, i.e. the lower bound is saturated. When applied to fracton models Chamon (2005); Bravyi et al. (2011); Haah (2011); Yoshida (2013); Vijay et al. (2015, 2016), the lower bound depends on the sizes and relative locations of the subsystems. For fracton models, there exist choices of subsystems for which the lower bound is nonzero and extensive. It is possible to have $S_{nonlocal}>0$ for subsystem choices, that are expected to have $S_{nonlocal}=0$ if TQFT holds.

This method observes an intimate relation between $S_{nonlocal}$ and topological excitations created in $D$ by a unitary operator $U$ stretched out in $CD$ which could be “deformed” into a unitary operator $U^{def}$ in $AD$, and since deformable $U$ is intimately related to condensate operator $W$ (for more details of $U$, $U^{def}$ and condensate operator $W$ see Sec.II.4), this method observes an intimate relation between $S_{nonlocal}$ and ground state condensates as well. This method allows a lower bound of $S_{nonlocal}$ to be obtained without calculating the entanglement entropy of any individual subsystem. Furthermore, it provides us with a unified viewpoint to understand the topology-dependent $S_{nonlocal}$ in the conventional topological orders and the geometry-dependent $S_{nonlocal}$ in fracton topological orders. Also discussed is the stability of the lower bound under local perturbations.

Two additional papers appeared after our work which also study the entanglement entropy of fracton phases, using explicit computation Ma et al. (2017b) and tensor network He et al. (2017).

For non-Abelian models, some of our assumptions breakdown, and our original method does not apply. Nevertheless, a variant of our lower bound is applicable to non-Abelian models Shi and Lu (2018).

The structure of the paper is as follows: In Sec.II we provide a derivation of the lower bound from some assumptions about topological excitations and condensates; In Sec.III we apply our lower bound to several exactly solved Abelian stabilizer models of 2D, 3D conventional topological orders and type I, type II fracton topological orders. In Sec.IV we discuss the stability of the lower bound under local perturbations. Sec.V is discussion and outlook.

II The lower bound

II.1 A few notations and definitions

We will consider a infinite system without boundaries. The system is divided into subsystems $A$, $B$, $C$, $D$ (nonoverlapping regions in real space, the union of which is the whole system). Each subsystem has a size large compared to the correlation length $\xi$, and the subsystems $A$ and $C$ are separated by a distance much larger than the correlation length. One example is shown in Fig.1, and similar constructions can apply to any dimensions. For all the examples in this paper, we have chosen $B$, $C$, $D$ to be local subsystems while $A$ is not, but there exist other possible choices, say $A$, $B$, $C$ local. A local subsystem is a subsystem which can be contained in a ball-shaped subsystem of finite radius. $\partial A$, $\partial B$ are the boundaries of the subsystems. We use $\bar{A}$ to denote the complement of $A$.

We will use $\rho$, $\sigma$ to represent density matrices. In this paper we always use $\rho$ for the ground state density matrix, $\rho=|\psi\rangle\langle\psi|$, and $|\psi\rangle$ is the ground state. We use $\rho_{ABC}$, $\sigma_{BC}$ when we want to specify the subsystems. The entanglement entropy is defined in terms of the (reduced) density matrix as usual $S=-\textrm{tr}[\sigma\ln\sigma]$. We use $S_{ABC}|_{\rho}$ and $S_{ABC}|_{\sigma}$ to distinguish the entanglement entropy on region $ABC$ with different density matrices $\rho$, $\sigma$. Define conditional mutual information

[TABLE]

and we use $I(A:C|B)|_{\rho}$ when we want to specify a density matrix. It is known that the conditional mutual information is always nonnegative $I(A:C|B)\geq 0$. We say $\sigma_{ABC}$ is conditionally independent if $I(A:C|B)|_{\sigma}=0$.

For unitary operators $U$ and $U^{\prime}$ which create excitations in $D$ when acting on the ground state $|\psi\rangle$, we say $U\sim U^{\prime}$ or $U$ is similar to $U^{\prime}$ if the states $U|\psi\rangle$ and $U^{\prime}|\psi\rangle$ have identical reduced density matrices on $ABC$, i.e. $tr_{D}[U\rho U^{\dagger}]=tr_{D}[U^{\prime}\rho U^{\prime\dagger}]$. Otherwise, we say $U$ and $U^{\prime}$ are distinct.

II.2 Prepare for the lower bound

If there is a density matrix $\sigma$ which is related to the ground state density matrix $\rho$ by $\rho_{AB}=\sigma_{AB}$ and $\rho_{BC}=\sigma_{BC}$, then we have:

[TABLE]

and the “$=$” happens if and only if $I(A:C|B)|_{\sigma}=0$. For a proof, observe that $I(A:C|B)|_{\rho}$ and $I(A:C|B)|_{\sigma}$ has only a single different term, and that $I(A:C|B)|_{\sigma}\geq 0$.

For a pure state, the entanglement entropy of a subsystem equals the entanglement entropy of its complement, e.g. $S_{\Omega}=S_{\bar{\Omega}}$ for any subsystem $\Omega$. Therefore:

[TABLE]

Observe that the local contributions of the entanglement entropy get canceled due to the fact that $A$ and $C$ are separated. Let us define the nonlocal entanglement entropy (of the ground state)

[TABLE]

The nonlocal entanglement entropy is just another way to write down the conditional mutual information, $S_{nonlocal}=I(A:C|B)|_{\rho}$, and therefore $S_{nonlocal}\geq 0$. The form in Eq.(2) has the advantage that it involves only local systems $B$, $C$, $D$. When the system is placed on a torus or other nontrivial manifolds instead of a infinite manifold, the system may have several locally indistinguishable ground states, this form of $S_{nonlocal}$ in terms of local subsystems is more convenient, and even if $\rho$ is a mixed state density matrix of different locally indistinguishable ground states, $S_{nonlocal}$ still has the same value.

II.3 The key idea about the lower bound

The discussion above suggests a way to obtain a lower bound of $S_{nonlocal}$. For any $\sigma$ satisfying $\sigma_{AB}=\rho_{AB}$ and $\sigma_{BC}=\rho_{BC}$,

[TABLE]

The density matrix $\sigma$ does not have to be a density matrix of a pure state. If we could find a $\sigma$ satisfying the above requirement and $S_{ABC}|_{\sigma}>S_{ABC}|_{\rho}$, a nonzero lower bound is obtained, and then the existence of nonzero nonlocal entanglement entropy is established.

Now, let us assume that we could find a set of $\sigma_{I}$ with $I=1,\ldots,N$ such that $\sigma_{I\,AB}=\rho_{AB}$ and $\sigma_{I\,BC}=\rho_{BC}$. Then we can do superpositions and define ${\sigma}\equiv\sum_{I=1}^{N}p_{I}\sigma_{I}$ with $\{p_{I}\}$ being a probability distribution, i.e. $p_{I}\in[0,1]$ and $\sum_{I=1}^{N}p_{I}=1$. The ${\sigma}$ is a new density matrix which satisfies ${\sigma}_{AB}=\rho_{AB}$ and ${\sigma}_{BC}=\rho_{BC}$. So we have a whole parameter space of ${\sigma}$ to try.

If lucky, we may even be able to find a $\sigma^{\ast}_{ABC}$ which is conditionally independent (satisfying $I(A:C|B)|_{\sigma^{\ast}}=0$) and we have an exact result

[TABLE]

Or, if we find the lower bound is saturated, we know the $\sigma$ we used to obtain the lower bound is conditionally independent. We note that, in the quantum case (unlike the classical case), it is not always possible to find a conditionally independent $\sigma^{\ast}$ such that $\sigma^{\ast}_{AB}=\rho_{AB}$ and $\sigma^{\ast}_{BC}=\rho_{BC}$ for $\rho$ being a general density matrix Ibinson et al. (2008). Therefore, the existence of such $\sigma^{\ast}$ in some system might be interesting by itself. On the other hand, the conditional independent state $\sigma^{\ast}$ is known to exist for models satisfying simple conditions (I)(II) in Kato et al. (2016).

II.4 Calculate the lower bound for Abelian models employing assumptions

We make a few assumptions about the condensates and operators creating topological excitations in order to develop a way to find $\sigma_{I}$ and calculate a lower bound of $S_{nonlocal}$. These assumptions are applicable to Abelian models with commuting projector Hamiltonians (with each term acting on a few sites localized in real space), including conventional models e.g. the toric code model in various dimensions Dennis et al. (2002), quantum double models with any Abelian finite group, and fracton models Chamon (2005); Bravyi et al. (2011); Haah (2011); Yoshida (2013); Vijay et al. (2015, 2016) of type I and type II which we will be focusing on in this work. In the context of fracton models, Abelian means no protected degeneracy associated with excitations.

Our way of employing the operators is inspired by a method by Kim and Brown Kim and Brown (2015), where an interesting connection between conditional mutual information and deformable operator $U$ is obtained. While the subsystems we choose have only an unimportant difference from Kim and Brown (2015), our result is different. The result in Kim and Brown (2015) shows that if $S_{nonlocal}=0$, there will be no topological excitations, and therefore, $S_{nonlocal}>0$ is needed for the existence of topological excitations. The result is very general since very small amount of assumptions was used. Nevertheless, the result was not powerful as a lower bound for $S_{nonlocal}$. In this work, on the other hand, we use more detailed properties of topological excitations and condensates to obtain a powerful lower bound of $S_{nonlocal}$. Our method shows that the key to have $S_{nonlocal}>0$ in these models is the nonlocal nature of the ground state condensates and the operators creating topological excitations. $S_{nonlocal}$ can be extensive in the subsystem size and it is not necessarily topologically invariant. Whether $S_{nonlocal}$ is topologically invariant or not depends on whether the operators can be deformed topologically.

Before rigorously stating the assumptions and deriving the lower bound, here are a few words about the physical picture. The ground states of topological orders condense extended objects. If a unitary operator $W$ (which has an extended support) acting on the ground state $|\psi\rangle$ gives you $W|\psi\rangle=e^{i\varphi}|\psi\rangle$, we call the operator $W$ a condensate operator with eigenvalue $e^{i\varphi}$. Whenever confusion can be avoided, we may call $W$ a condensate for short. Let us further assume $W$ to be a tensor product of operators acting on each site. Then, a suitably defined “truncation” of a condensate operator $W$ onto a subsystem $\Omega$ gives you a new operator $U$. $U|\psi\rangle$ is an excited state with topologically excitations located around $\partial\Omega\cap W$, and $U$ can be deformed in the sense that you could choose a “truncation” of $W^{\dagger}$ onto $\bar{\Omega}$ and call it $U^{def}$, which creates the same topological excitations and satisfies $U|\psi\rangle=U^{def}|\psi\rangle$. One the other hand, if we have unitary operators $U$ and $U^{def}$ satisfying $U|\psi\rangle=U^{def}|\psi\rangle$, then $[U^{def}]^{\dagger}U|\psi\rangle=|\psi\rangle$, and therefore $[U^{def}]^{\dagger}U$ is a condensate operator with eigenvalue $+1$. Intuitively, a condensate operator $W$ can be “truncated” into a deformable operator $U$ which creates topological excitations, and a pair of deformable operators $U$ and $U^{def}$ can be “glued together” into a condensate operator $W$. Therefore $U$ and $W$ are closely related. Some examples of condensate operator $W$ and deformable operator $U$, are shown in Fig.2. Once we have condensate operators $W_{i}$ and deformable operators $U_{j}$, we use $U_{j}$ to create topological excitations in $D$ which result in some states $\sigma_{I}$ which is identical to the ground state $\rho$ on $AB$ and $BC$. If the excitations created in $D$ are topological excitations, they can not be created by an operator supported on $D$, and $\sigma_{I\,ABC}$ and $\rho_{ABC}$ will have some difference. The difference is detected by a change of the eigenvalue of condensate operator $W_{i}$ supported on $ABC$. Then, we use $\sigma_{I}$ to obtain a lower bound of $S_{nonlocal}$.

The following are our assumptions $\mathbf{U}$-1, $\mathbf{U}$-2, $\mathbf{U}$-3; $\mathbf{W}$-1, $\mathbf{W}$-2:

Assumption $\mathbf{U}$-1: There exists a set of unitary operators $\{U_{i}\}$ supported on $CD$ and a set of unitary operators $\{{U}_{i}^{def}\}$ supported on $AD$, with $i=1,\cdots,M$. See Fig.3 for an example. When acting on the ground state $|\psi\rangle$, $U_{i}$ can be “deformed” into $U_{i}^{def}$, i.e.:

[TABLE]

Assumption $\mathbf{U}$-2: For any subsystem $\Omega$, the unitary operator $U_{i}$ can always be written as a direct product of unitary operators $U_{i\,\Omega}$ and $U_{i\,\bar{\Omega}}$ which act on the subsystem $\Omega$, $\bar{\Omega}$ respectively:

[TABLE]

Assumption $\mathbf{U}$-3: There are integers $n_{i}$ such that $U_{i}^{n_{i}}=1$, and when multiple $U_{i}$ act on a ground state $|\psi\rangle$, we have the following:

[TABLE]

where integer $0\leq k_{i}\leq n_{i}-1$, and we allow possible phase factors $e^{i\delta(\{k_{i}\})}$.

Assumption $\mathbf{W}$-1: There exists a set of unitary operators {$W_{i}$} supported on subsystem $ABC$ such that

[TABLE]

Assumption $\mathbf{W}$-2: The following relation between $W_{i}$ and $U_{j}$ holds:

[TABLE]

where $\delta_{ij}$ is the Kronecker delta.

Comments about the assumptions:

1. $\mathbf{U}$-1 implies that when $U_{i}$ is acting on the ground state, it can create excitations only in $D$, but not in $ABC$.

2. We do not assume $U_{i}$, $W_{i}$ to be string operators and not even assume $U_{i}$, $W_{i}$ to be integer dimensional operators. In fact, we will apply this method to fractal operators later.

3. In $\mathbf{U}$-1 we assumed $U$ can be deformed without changing the excitations. Nevertheless, we do not assume $U$ can be deformed into all topologically equivalent configurations, and we do not assume $U$ can be deformed continuously. As we will see below, in fracton models, deformations exist in weird form, $U$ may not be topologically deformed (i.e. deformed continuously into any topologically equivalent configuration), and $U$ can sometimes be deformed in discontinuous ways into topologically inequivalent configurations.

4. $\mathbf{U}$-2 is not true when a local perturbation is added. We will address the stability of the lower bound under local perturbations separately in Sec.IV.

5. For the non-Abelian case, the entanglement entropy of an excited state could depend on the quantum dimension of the anyon Kitaev and Preskill (2006); Shi and Lu (2018), and $\mathbf{U}$-2 does not apply. On the other hand, the idea in Sec.II.3 still holds, a saturated lower bound for non-Abelian models is recently discussed in Shi and Lu (2018).

6. For systems with boundaries, one may choose $D$ being region attached to boundaries, as is done in Kim and Brown (2015). An alternative way is to identify $D$ with a boundary region $\partial\Omega$; in this case, $\mathbf{U}$-1 should be understood as: $U_{i}$ being an operator supported on $C$ and attached to $\partial C\cap\partial\Omega$, and $U_{i}^{def}$ being an operator supported on $A$ and attached to $\partial A\cap\partial\Omega$.

7. According to $\mathbf{W}$-2, $|\psi\rangle$ and $U_{i}|\psi\rangle$ are eigenstates of $W_{i}$ with different eigenvalues, where $|\psi\rangle$ is the ground state. Since $W_{i}$ is supported on $ABC$, this implies that $U_{i}$ is distinct from the identity operator. Similarly, $U_{i}$ and $U_{j}$ are distinct for $i\neq j$. We will refer to this change of eigenvalue of $W_{i}$ as a detection, e.g. $U_{i}$ is detected by $W_{i}$. The requirement $\theta_{ij}=\frac{2\pi}{n_{i}}\delta_{ij}$ is not crucial, and it can be replaced by other numbers as long as the operator set $\{W_{i}\}$ can detect the difference among the set of operators $\{U_{i}\}$.

8. For a relatively simple class of models, which has the Hamiltonian $H=-\sum_{k}h_{k}$, $[h_{i},h_{j}]=0$ and $h_{k}^{2}=1$, there is an obvious class of operators that satisfy Eq.(8) in $\mathbf{W}$-1, namely $W_{i}=\prod_{k\in\mathcal{E}_{i}}h_{k}$, where $\mathcal{E}_{i}$ is a subset of the stabilizer generators. For the ground state $|\psi\rangle$, we have $h_{i}|\psi\rangle=|\psi\rangle$, it follows that $W_{i}|\psi\rangle=|\psi\rangle$. As $U_{i}$ could not flip stabilizers in $ABC$, $\mathcal{E}_{i}$ must contain some $h_{j}$ in $D$. It turns out that this simple observation applies to all the Abelian stabilizer models we will use as examples in Sec.III. However, we do not provide a general procedure to find the subset $\mathcal{E}_{i}$ for fracton models. On the other hand, our method works for models not in this simple class also, such as quantum double models Kitaev (2003); Bombin and Martin-Delgado (2008) with Abelian finite groups.

Define the set of states

[TABLE]

with $k_{i}$ being integers. Define $\sigma(\{k_{i}\})\equiv|\{k_{i}\};\psi\rangle\langle\{k_{i}\};\psi|$. Note the total number of $\sigma(\{k_{i}\})$ is $\prod_{i=1}^{M}n_{i}$. Relabel $\sigma(\{k_{i}\})$ using a new index $I=1,\cdots,N$ with $N=\prod_{i=1}^{M}n_{i}$ and call them $\sigma_{I}$. One immediately varifies that

1. $\mathbf{U}$-1, $\mathbf{U}$-3 $\Rightarrow$ $\sigma_{I\,AB}=\rho_{AB}$ and $\sigma_{I\,BC}=\rho_{BC}$;

2. $\mathbf{W}$-1, $\mathbf{W}$-2 $\Rightarrow$ $\sigma_{I\,ABC}\cdot\sigma_{J\,ABC}=0$ for $I\neq J$;

3. $\mathbf{U}$-2 $\Rightarrow$ $\sigma_{I\,ABC}=V_{I}\rho_{ABC}V_{I}^{\dagger}$ and $S_{ABC}|_{\sigma_{I}}=S_{ABC}|_{\rho}$.

Where $V_{I}$ is some unitary operator acting on subsystem $ABC$ and recall that $\rho$ is the ground state density matrix.

Let ${\sigma}=\sum_{I=1}^{N}p_{I}\sigma_{I}$ with probability distribution $\{p_{I}\}$, one derives that

[TABLE]

“=” if and only if $p_{I}=\frac{1}{N}$ for all $I$. From Eq.(3) we find

[TABLE]

II.5 When is the lower bound topologically invariant?

It is instructive to think of the conditions under which our lower bound of $S_{nonlocal}$ is topologically invariant.

Consider a chosen set of subsystems $A$, $B$, $C$, $D$ and the operator sets $\{U_{i}\}$, $\{U_{i}^{def}\}$ and $\{W_{i}\}$. Let us do “topological deformations” of the subsystems and the operator sets. Here, by “topological deformation” of the operator sets we mean that we can topologically deform the support of each operator to get new operator sets which preserve the algebra in Eq.(5,6,7,8,9). Note that, these deformations generally change the positions of the excitations, which should be contrasted with the type of deformation in $\mathbf{U}$-1, in which the positions of the excitations never change. When these conditions are satisfied, the lower bound for the two topologically equivalent choices of subsystems are the same. If such conditions are satisfied for each pair of topologically equivalent choices of subsystems, then our lower bound will be topologically invariant.

As is shown in the examples below in Sec.III, $S_{nonlocal}$ can be either topologically invariant or not, and it is instructive to think of how the conditions above are violated in fracton models Chamon (2005); Bravyi et al. (2011); Haah (2011); Yoshida (2013); Vijay et al. (2015, 2016) in which $S_{nonlocal}$ depends on the geometry of subsystems.

III Applications

In this section, our lower bound is applied to several stabilizer models of Abelian phases: the 2D and 3D conventional topological orders and type I, type II fracton phases.

III.1 The 2D Toric Code Model

For a 2D topological order, choose the subsystems $A$, $B$, $C$, $D$ of the same topology as is shown in Fig.1. From the well-known results Kitaev and Preskill (2006); Levin and Wen (2006), one derives $S_{nonlocal}=2\gamma$ where $-\gamma$ is the topological entanglement entropy.

For the 2D toric code model Kitaev (2003): On a square lattice with a qubit on each link, the Hamiltonian is

[TABLE]

where $A_{s}$ is a product of $X_{r}$ of a “star” or vertex, $B_{p}$ is a product of $Z_{r}$ of a plaquette.

[TABLE]

where $X_{r}$, $Z_{r}$ are Pauli operators acting on the qubit on link $r$.

The ground state of toric code model condenses two types of closed string operators, and the corresponding open string operators (which could be regarded as truncations of closed string operators) create topological excitations at the endpoints.

We find the following unitary operators $U_{1}$, $U_{2}$, $W_{1}$, $W_{2}$ as is shown in Fig.4. $U_{1}$ and $W_{2}$ are products of $X_{r}$; $U_{2}$ and $W_{1}$ are products of $Z_{r}$. Also notice the feature that the closed string operators $W_{1}$ ($W_{2}$) can be written as a product of stabilizers $B_{p}$ ($A_{s}$) on a 2D disk region surrounded by the corresponding closed strings.

$\mathbf{U}$-1, $\mathbf{U}$-2, $\mathbf{U}$-3, $\mathbf{W}$-1, $\mathbf{W}$-2 can be checked. The operators satisfy:

[TABLE]

Therefore $M=2$, $n_{1}=n_{2}=2$ and $N=n_{1}n_{2}=4$. Using the result in Eq.(10), one derives $S_{nonlocal}\geq 2\ln 2$. By comparing with the known result $\gamma=\ln 2$, $S_{nonlocal}=2\gamma=2\ln 2$, we find that our lower bound is saturated.

A by-product of a saturated lower bound is an explicit construction of a conditionally independent $\sigma^{\ast}$. In the toric code case:

[TABLE]

and $\sigma^{\ast}$ satisfies:

1. $\sigma^{\ast}_{AB}=\rho_{AB}$, $\sigma^{\ast}_{BC}=\rho_{BC}$;

2. $I(A:C|B)|_{\sigma^{\ast}}=0$.

Where $\rho=|\psi\rangle\langle\psi|$ is the ground state density matrix.

The observation in Sec.II.5 explains why the lower bound is topologically invariant in the toric code model: the operators $U_{i}$ and $W_{i}$ can be topologically deformed together with the subsystems $A$, $B$, $C$ and $D$, without changing the algebra in Eq.(5,6,7,8,9).

This method can be applied to other 2D Abelian topological orders, e.g., quantum double models with Abelian finite groups, and the lower bounds are saturated. For a variant of the method for non-Abelian models, see Shi and Lu (2018).

III.2 The 3D Toric Code model

The 3D toric code model Hamma et al. (2005a) is defined on a cubic lattice, with one qubit on each link. The Hamiltonian is of exactly the same form as the one of the 2D toric code model:

[TABLE]

[TABLE]

Here a star $s$ includes the 6 links around a vertex, and a plaquette $p$ is a square consistent of 4 links.

The ground state of the 3D toric code model condenses one type of closed string and one type of closed membrane. There is one type of open string operator that creates point-like topological excitations at the endpoints and one type of open membrane operator which creates loop-like topological excitations at the edge of the membrane, see Fig.5.

III.2.1 Subsystem types for 3D models

The 3D toric code model is the first 3D model we discuss, and it is a good place to introduce subsystem types for 3D models which will be discussed for all 3D models. We focus on the following three topologically distinct subsystem types, i.e. the type-$\alpha$,$\beta$,$\gamma$ shown in Fig.6, although other choices are possible. We will use the notation $S^{(\alpha)}_{nonlocal}$, $S^{(\beta)}_{nonlocal}$, $S^{(\gamma)}_{nonlocal}$ to distinguish the nonlocal entanglement entropy for the three topological types.

Type-$\alpha$ has $D$, which consists of two disconnected boxes, while $CD$ is connected, and it can be used to detect open string-like $U_{i}$, which is attached to the two boxes of $D$. Type-$\beta$ has $D$ of the topology of a solid torus (and therefore $D$ is not simply connected), while $CD$ is simply connected. It can be used to detect open membrane-like $U_{i}$ supported on $CD$ which create excitations in $D$ as noncontractible loops. Type-$\gamma$ has $C$ and $CD$ of the same topology, e.g. the topology of a solid torus, and $B$ is simply connected.

Type-$\alpha$ and type-$\beta$ have already been implied in paper by Kim and Brown Kim and Brown (2015), in which, similar subsystems types are used to study different types of boundaries of 3D models.

For type-$\gamma$, $S_{nonlocal}^{(\gamma)}=0$ for models satisfying the assumptions in Grover et al. (2011), e.g. the entanglement entropy of a general subsystem $\Omega$ (which has large size compared to correlation length) can be decomposed into local plus topological parts:

[TABLE]

Therefore, a model with $S_{nonlocal}^{(\gamma)}>0$ is a model beyond the description of Grover et al. (2011). Fractal models do have $S^{(\gamma)}_{nonlocal}>0$ for some choices of the subsystems and the value can be extensive. Furthermore, the contribution of $S^{(\alpha)}_{nonlocal}$ ($S^{(\beta)}_{nonlocal}$) is not necessarily from open string-like (open membrane-like) $U$.

III.2.2 The 3D Toric Code model has saturated lower bounds for each subsystem type

Let us go back to the 3D toric code model.

For type-$\alpha$, we find $M=1$, $n_{1}=2$ where the operator $U_{1}$ is an open string operator which creates a point-like topological excitation in each box of $D$ and $W_{1}$ is a closed membrane operator. Therefore $N=2$ and $S_{nonlocal}^{(\alpha)}\geq\ln 2$.

For type-$\beta$, we find $M=1$, $n_{1}=2$ where the operator $U_{1}$ is an open membrane operator which creates a loop-like topological excitation at the edge of the membrane (the loop could not continuously shrink within $D$ into a point), and $W_{1}$ is a closed string operator. Therefore $N=2$ and $S_{nonlocal}^{(\beta)}\geq\ln 2$.

For type-$\gamma$, operators supported on $CD$, which create excitations in $D$, could always be deformed into $D$. This is because, for 3D toric code, the operators that create topological excitations can be topologically deformed keeping the excitations fixed. Therefore we obtain a lower bound $S_{nonlocal}^{(\gamma)}\geq 0$.

Comparing with the known entanglement properties Grover et al. (2011); Kim and Brown (2015) of the 3D toric code model, our lower bounds are identical to the exact results:

[TABLE]

III.3 The X-Cube Model

The X-cube model is a 3D exactly solved stabilizer model, and it is an example of type I fracton phase Vijay et al. (2016). The model is defined on a cubic lattice with one qubit on each link. The Hamiltonian is

[TABLE]

where $A_{c}$ is a product of $Z_{r}$ on a cube (which includes 12 links), and $B^{(xy)}_{s}$, $B^{(yz)}_{s}$, and $B^{(zx)}_{s}$ are products of $X_{r}$ of 4 links around a vertex which are parallel to the $xy$-plane, $yz$-plane, $zx$-plane respectively.

Here we focus on type-$\alpha$ and consider the geometry dependence of $S_{nonlocal}$, see Fig.7. We find the following lower bounds:

Fig.7a, $S^{(\alpha)}_{nonlocal}\geq(2l_{D}+O(1))\ln 2$;

Fig.7b, $S^{(\alpha)}_{nonlocal}\geq(4l_{D}+O(1))\ln 2$;

Fig.7c, $S^{(\alpha)}_{nonlocal}\geq 0$.

Where $O(1)$ denotes order one contributions which dependent on the detailed shapes of the subsystems, which is not crucial for our discussion.

The types of $U_{i}$ that contribute to the $S^{(\alpha)}_{nonlocal}$ in Fig.7a are illustrated in Fig.8a, each $U_{i}$ stretches out in directions parallel to the $yz$-plane and creates a pair of “dimension-2 anyons.” The translations of the operators $U_{i}$ in Fig.8a in $(1,0,0)$ direction give you distinct operators (while translations in $(0,1,0)$ or $(0,0,1)$ directions do not give you distinct operators). This gives $M=2l_{D}+O(1)$.

Translations can produce distinct operators, this indicates a breakdown of topological deformation: in the X-cube model, $U_{i}$ is deformable but not topologically deformable. Another nice example of the breakdown of topological deformation is shown in Fig.8c, in which $U$ and its translations $U_{T_{1}}$, $U_{T_{2}}$, and $U_{T_{1}+T_{2}}$ are distinct; nevertheless, by thinking of the condensate in Fig.8b, one can show $U\sim U_{T_{1}}\cdot U_{T_{2}}\cdot U_{T_{1}+T_{2}}$.

The result $S^{(\alpha)}_{nonlocal}\geq(4l_{D}+O(1))\ln 2$ for Fig.7b can be understood by thinking of contributions from operators parallel to the $yz$-plane and the $xy$-plane, which gives $M=4l_{D}+O(1)$. The result $S^{(\alpha)}_{nonlocal}\geq 0$ in Fig.7c comes from the fact that the types of $U_{i}$ discussed above could not connect the two boxes of $D$ separated by a displacement vector $\vec{d}=(d_{x},d_{y},d_{z})$ with $|d_{x}|,|d_{y}|,|d_{z}|>l_{D}$ and the inability to find $U_{i}$ gives $M=0$.

These results for $S^{(\alpha)}_{nonlocal}$ may also be calculated using the method in Hamma et al. (2005b) and an independent estimation agrees with our lower bounds up to $O(1)$ contributions.

The lower bounds of $S^{(\alpha)}_{nonlocal}$ for all the cases in Fig.7 depend on the length scale $l_{D}$ and the displacement vector $\vec{d}$ but are not sensitive to other details. This is due to the fact that we have chosen “big enough” $A,B,C$, so that they do not block any $U_{i}$. If we consider another extreme, say the subsystem $C$ has a very narrow neck, then $S^{(\alpha)}_{nonlocal}$ will be sensitive to the geometry of the neck which determines how many $U_{i}$ could pass through.

One may also apply the same idea to subsystems of type-$\beta$ and type-$\gamma$ and find extensive values of $S^{(\beta)}_{nonlocal}$ and $S^{(\gamma)}_{nonlocal}$ for certain choices of subsystems.

III.4 Fractal spin liquids

Fractal spin liquids Yoshida (2013) is a generalization of Haah’s code Haah (2011). A common feature of fractal spin liquid models is the existence of fractal condensates. Fractal structures have discrete scale symmetries, and this results in a more complicated dependence of the ground state degeneracy on the system size Haah (2013); Yoshida (2013) compared to the type I fracton models.

Some of the fractal models possess “hybrid” condensates $W_{i}$ having both 1D parts and fractal parts, and the truncations of the condensates give you $U_{i}$ which can be either a string-like operator or a fractal operator. Note that, they do not fit into the definition of type I due to the existence of fractal operators. On the other hand, they do not fit into type II because the excitations created by the string-like operator are mobile excitations. Some fractal models have only fractal condensates, and no string-like $U_{i}$ exists. These models are type II fracton models.

Discussed in the following are ways to detect string-like $U_{i}$ and fractal $U_{i}$ using condensates. Then, $S_{nonlocal}$ is shown to be extensive for certain choices of subsystem geometry.

III.4.1 The Sierpinski Prism Model

As an example, we consider the model (d) in Yoshida’s paper Yoshida (2013). Let us call this model the Sierpinski prism model, named after the shape of the condensate in Fig.9a, which looks like a prism with three legs decorated with Sierpinski triangles. This model lives on a 3D cubic lattice with two qubits ($A$ and $B$) on each site. The Hamiltonian can be written as

[TABLE]

where $i,j,k$ are integers labeling the sites on cubic lattice, and the Hamiltonian involves all the translations of the operator $h^{Z}_{(0,0,0)}$ and $h^{X}_{(0,0,0)}$. Explicitly, in terms of Pauli operators acting on each $A$ and $B$ qubit on different sites, we have

[TABLE]

It is easy to check that all terms in the Hamiltonian commute and $[h_{(i,j,k)}^{Z}]^{2}=[h^{X}_{(i,j,k)}]^{2}=1$.

Using Yoshida’s notation:

[TABLE]

Where $\bar{x}\equiv x^{-1}$, $\bar{y}\equiv y^{-1}$ and $\bar{z}\equiv z^{-1}$. In terms of polynomials $f(x)$, $g(x)$ with coefficients over $\mathbf{F}_{2}$, i.e. the coefficients can take 0 or 1:

[TABLE]

Here $f(x)=1+x$ and $g(x)=1$ for the Sierpinski prism model. The polynomials with $\mathbf{F}_{2}$ coefficients indicate the locations and the numbers of Pauli $Z$ or $X$ operators in the product; the upper row is for $A$ qubits and the lower row is for $B$ qubits.

Choosing other polynomials $f(x)$, $g(x)$ or changing $\mathbf{F}_{2}$ into $\mathbf{F}_{p}$ ($p>2$ prime number) will generally give you other fractal models.

The Sierpinski prism model possesses hybrid condensates which consist of 1D parts and fractal parts, see Fig.9. The condensates in Fig.9a and Fig.9b can be constructed as a product of $h^{Z}_{(i,j,k)}$ and the condensate in Fig.9c can be constructed as a product of $h^{X}_{(i,j,k)}$. As is suggested by the discrete scaling symmetry of fractal structure and the continuous scaling symmetry of a 1D line: the upper and lower surfaces can be separated by an arbitrary distance in $z$-direction (without changing the size of upper/lower surfaces), and under a rescaling $l\to 2^{m}l$ the condensates look similar. Under other rescaling factors, the condensates look different but they could be constructed using a product of condensates that looks similar to the ones in Fig.9.

While this model does not have any logical qubits under periodic boundary conditions on an $L_{x}\times L_{y}\times L_{z}$ lattice, i.e. $x^{L_{x}}=y^{L_{y}}=z^{L_{z}}=1$, it does have logical qubits under some “twisted” boundary conditions (say $x^{L_{x}}=1$, $y^{L_{y}}=x$, $z^{L_{z}}=1$ with $L_{x}=2^{m}+1$, $L_{y}=2^{m}$, $\forall\,L_{z}$ and integer $m$), or under open boundary conditions.

Despite the fact that the Sierpinski prism model is one of the simplest fractal models, it nicely illustrates all the important ingredients needed in order to understand how our method works in fractal models. To be specific, it illustrates the following three types of detections:

1. The detection of a string-like $U$ using a fractal $W$.

2. The detection of a fractal $U$ using a string-like $W$.

3. The detection of a fractal $U$ using a fractal $W$.

III.4.2 The detection of sting-like $U$ by fractal $W$

As is shown in Fig.9, we have condensates with string-like parts and fractal parts. String-like $U_{i}$ can be obtained from a truncation of the condensates. In Fig.10a, translations of a string-like $U$, e.g. $U_{T_{1}}$ and $U_{T_{2}}$ (with $\vec{T}_{1}=(-2^{m},0,0)$, $\vec{T}_{2}=(0,2^{m},0)$ and an integer $m$) are distinct from $U$, but $U\sim U_{T_{1}}\cdot U_{T_{2}}$. It is different from what happens in conventional topological orders but similar to what happens in type I fracton models, see Fig.8.

The distinctness of a string-like $U$ and its translations indicates a lower bound of $S_{nonlocal}$ extensive in the subsystem size, and indeed we can use the fractal part of $W_{i}$ to detect different string-like $U_{i}$, see Fig.10b, and get an extensive lower bound.

III.4.3 The detection of fractal $U$ by string-like $W$

Very similarly, fractal $U$ can be detected by string-like parts of $W$, see Fig.11a. It is clear that a translation of a fractal $U$ will be a distinct operator if it anticommutes with a different $W$. Therefore, by suitably choosing the geometrical shapes of subsystems $A,B,C,D$, it is possible to get an extensive lower bound of $S_{nonlocal}^{(\alpha)}$. $S_{nonlocal}^{(\beta)}$ and $S_{nonlocal}^{(\gamma)}$.

III.4.4 The detection of fractal $U$ by fractal $W$

Another way to detect fractal $U$ is to use a fractal part of a condensate $W$, see Fig.11b, and it is the only way to detect $U$ for those fractal models without string operators i.e., when the condensates contain only a fractal structure without string-like parts (i.e., type II). Therefore, it is important to understand this case.

The key features, which can be observed in Fig.11b are the following.

1. The fractal condensates are supported on 2D surfaces with “holes” of different (discrete) length scales, in other words it is less than 2D.

2. Fractal $U$ and fractal (part of) $W$ lie in distinct intersecting surfaces, and it is possible to make the operators intersecting at a point. Certain translation of $U$ has a non-overlapping support with $W$ and therefore it commutes with $W$. This implies that translations of $U$ can give you distinct operators.

After some thought, one finds it is possible to get an extensive lower bound of $S_{nonlocal}^{(\alpha)}$, $S_{nonlocal}^{(\beta)}$ and $S_{nonlocal}^{(\gamma)}$ by suitably choosing the geometrical shapes of the subsystems $A,B,C,D$. A choice of $D$ for type-$\beta$ or type-$\gamma$ is shown in yellow color in Fig.11b.

III.4.5 Further comments

1. One may consider other subsystem types, for example: type-$\delta$ with $D$ consists of three disconnected boxes, and suitably chosen $A,B,C$. When putting the three boxes on the positions of the three excitations of a $U$ in Fig.11, $S^{(\delta)}_{nonlocal}>0$.

2. When rescaling the subsystem sizes according to the discrete scale symmetries (for the Sierpinski prism model, it is $L\rightarrow 2^{m}L$), the change of $S_{nonlocal}$ could be investigated using entanglement renormalization group transformation Haah (2014). For a rescaling by a factor which is not in the discrete scaling group, $S_{nonlocal}$ may change in more complicated way.

3. It is possible to have a fractal model with a unique ground state on a $T^{3}$, in which case, there still exists nonzero $S_{nonlocal}$. This indicates that $S_{nonlocal}$ is in some sense more universal than the ground state degeneracy.

4. For models with only fractal condensates (i.e. type II), although it is possible to find $S_{nonlocal}^{(\alpha)}>0$ for some choices of $A,B,C,D$, we get a lower bound [math] when, for example, the two boxes of $D$ of length $l_{D}\times l_{D}\times l_{D}$ are separated by a displacement vector $\vec{d}$ satisfying $|\vec{d}|>\lambda l_{D}$ with some constant $\lambda$ depending on the model. It might be a general exact result that $S^{(\alpha)}_{nonlocal}=0$ when $|\vec{d}|>\lambda l_{D}$ no matter how you choose $A,B,C$ but we do not have a proof. The case for Haah’s code is a conjecture by Kim 111I. H. Kim, private communication. . If the conjectured results are true, this may be used to make a clear distinction between type I and type II fracton models.

IV Perturbations

The stability of quantities under local perturbations is an extremely important topic. If some property of an exactly solved model is totally changed when a tiny local perturbation is added, this property could never be observed in real systems.

The ground state degeneracy of topological orders is robust (stable) to arbitrary local perturbations. It is known that the toric code model is stable under arbitrary local perturbations Kitaev (2003). The stability of ground state degeneracy is proved Bravyi et al. (2010); Bravyi and Hastings (2011) for a very general class of models which satisfy assumptions TQO-1 and TQO-2. In the proof, Osborne’s modification Osborne (2007) of the quasi-adiabatic continuation Hastings and Wen (2005) is employed. This proof is applicable to both conventional and fracton topological orders.

We would like to understand the stability of $S_{nonlocal}$ under local perturbations. It turns out that the stability of $S_{nonlocal}$ is a trickier problem compared to the stability of the ground state degeneracy. The corresponding problem for the conventional topological orders, e.g. the stability of topological entanglement entropy is not solved completely without additional assumptions. It is known from Bravyi’s counterexample (see Zou and Haah (2016) for a published reference) that the arguments provided in the original works Kitaev and Preskill (2006); Levin and Wen (2006) about the invariance of the topological entanglement entropy under perturbation are not complete. Kim obtained a bound of the change of topological entanglement entropy with a $1$st order perturbation Kim (2012) assuming the conditionally independence of certain subsystems.

Here we study the stability of our lower bound of $S_{nonlocal}$ under a finite depth quantum circuit for simplicity, since it is known from the viewpoint of quasi-adiabatic evolution Hastings and Wen (2005) that local perturbations for gapped systems can be approximated by finite depth quantum circuit Chen et al. (2010); Haah (2016).

Assumption $\mathbf{S}$: For subsystems $A$, $B$, $C$ as is shown in Fig.12. $\delta Q$ is a unitary operator which has support intersecting with $A$ and $C$. $B$ is separated from $\delta Q$ by a distance $d$. $\sigma$ is a density matrix of a state with correlation length $\xi$ and replica correlation length Zou and Haah (2016) $\xi_{\alpha}$ and $\sigma^{\prime}\equiv\delta Q\sigma\delta Q^{\dagger}$. And $\sigma$ is a density matrix such that

[TABLE]

where “$\simeq$” means there is a correction that is negligible when $d$ is large compared to $\xi$ and $\xi_{\alpha}$.

The assumption $\mathbf{S}$ should be understood as an assumption about the density matrix $\sigma$. It is trivial to check that “$\simeq$” can be replaced by “$=$” when $\sigma_{AB}=\sigma_{A}\otimes\sigma_{B}$, and it may seem intuitive that the difference between the left-hand side and the right-hand side of Eq.(14) should decay as $e^{-d/\xi}$. Nevertheless, the original suggestion Kitaev and Preskill (2006) that $\mathbf{S}$ is true for $\xi\ll d$ is violated in Bravyi’s counterexample. It is observed in Zou and Haah (2016) that this is due to the fact that the replica correlation length $\xi_{\alpha}$ is infinity for the cluster state in Bravyi’s counterexample. When the cluster state is deformed, $\xi_{\alpha}$ becomes finite. For generic local perturbations without symmetry requirement, it is fine-tuned to have $\xi_{\alpha}=+\infty$ but $\xi_{\alpha}$ can be arbitrarily large compared to $\xi$. Judging from a recent conjecture Zou and Haah (2016) , $\xi_{\alpha}\ll d$ may be the condition required for $\mathbf{S}$ to be true.

In the following, we discuss the stability of our lower bound of $S_{nonlocal}$ under a depth-$R$ quantum circuit $Q$ which creates a perturbed ground state $\bar{\rho}$ satisfying $\mathbf{S}$. We take $R\sim\xi$, and assume $\xi$ and $\xi_{\alpha}$ much smaller than the length scales of the subsystems. This analysis does not cover all possible local perturbations (especially those with $\xi_{\alpha}\to+\infty$), but we believe it covers a large class of interesting local perturbations.

Let $Q$ be the depth-$R$ quantum circuit ($QQ^{\dagger}=1$) which is responsible for the local perturbation. In other words, we assume the following objects in the perturbated model are related to the corresponding objects in the unperturbed model by

1. The new Hamiltonian: $\bar{H}=QHQ^{\dagger}$;

2. The new (dressed) operators: $\bar{U}_{i}=QU_{i}Q^{\dagger}$, $\bar{U}_{i}^{def}=Q\bar{U}_{i}^{def}Q^{\dagger}$, and $\bar{W}_{i}=QW_{i}Q^{\dagger}$;

3. The new ground state: $|\bar{\psi}\rangle=Q|\psi\rangle$;

4. The new density matrices: $\bar{\sigma}_{I}=Q\sigma_{I}Q^{\dagger}$ and $\bar{\rho}=Q\rho Q^{\dagger}$.

The dressed operators typically have a “fatter” support than the corresponding operators in the unperturbed stabilizer model, see Fig.(13) for an illustration. It is possible that the support of some $\bar{U}_{i}$ will overlap with $AB$, the support of some $\bar{U}_{i}^{def}$ will overlap with $BC$ and the support of some $W_{i}$ will overlap with $D$. For those operators, we need to throw them away, and supply with other operators if possible. We label the remaining operators using $\bar{i}$, $\bar{i}=1,\cdots\bar{M}$, where $\bar{M}\leq M$, and we call the remaining density matrices $\bar{\sigma}_{\bar{I}}$, with $\bar{I}=1,\cdots,\bar{N}$, where $\bar{N}=\prod_{\bar{i}=1}^{\bar{m}}\bar{n}_{\bar{i}}$.

With $\mathbf{U}$-1, $\mathbf{U}$-2, $\mathbf{U}$-3, $\mathbf{W}$-1, $\mathbf{W}$-2 satisfied for the unperturbed system, we supply with $\mathbf{S}$ in order to complete a result about the stability of the lower bound. With these assumptions, one could verify the following results:

$1^{\prime}$) $\mathbf{U}$-1, $\mathbf{U}$-3 $\Rightarrow$ $\bar{\sigma}_{\bar{I}\,AB}=\bar{\rho}_{AB}$ and $\bar{\sigma}_{\bar{I}\,BC}=\bar{\rho}_{BC}$;

$2^{\prime}$) $\mathbf{W}$-1, $\mathbf{W}$-2 $\Rightarrow$ $\bar{\sigma}_{\bar{I}\,ABC}\cdot\bar{\sigma}_{\bar{J}\,ABC}=0$ for $\bar{I}\neq\bar{J}$;

$3^{\prime}$) $\mathbf{U}$-1, $\mathbf{U}$-2, $\mathbf{U}$-$3$, $\mathbf{S}$ $\Rightarrow$ $S_{ABC}|_{\bar{\sigma}_{\bar{I}}}\simeq S_{ABC}|_{\bar{\rho}}$.

The derivation of $1^{\prime}$) and $2^{\prime}$) are parallel to what is done in Sec.II.4. The derivation of $3^{\prime}$) follows. There exists a unitary operator $\delta Q$ supported on a region within a distance $R$ around $\partial C\cap\partial D$, See Fig.14, such that $\bar{\rho}=\delta Q\bar{\rho}^{\prime}\delta Q^{\dagger}$, $\bar{\sigma}_{\bar{I}}=\delta Q\bar{\sigma}^{\prime}_{\bar{I}}\delta Q^{\dagger}$ and $\bar{\sigma}^{\prime}_{\bar{I}\,ABC}=\bar{V}_{\bar{I}}\bar{\rho}^{\prime}_{ABC}\bar{V}^{\dagger}_{\bar{I}}$. $\bar{V}_{\bar{I}}$ is some unitary operator supported on $ABC$. It is always possible to find such $\delta Q$ and $\bar{V}_{\bar{I}}$ given that $\mathbf{U}$-2 is satisfied for the unperturbed case. Therefore, $S_{ABC}|_{\bar{\sigma}^{\prime}_{\bar{I}}}=S_{ABC}|_{\bar{\rho}^{\prime}}$.

Then, by applying $\mathbf{U}$-$1$, $\mathbf{U}$-$3$, and $\mathbf{S}$ we find that

[TABLE]

After simple algebra one arrives at the result $3^{\prime}$) i.e. $S_{ABC}|_{\bar{\sigma}_{\bar{I}}}\simeq S_{ABC}|_{\bar{\rho}}$.

Assuming the error caused by “$\simeq$” could be neglected, one could apply the same method as Sec.II.4 to arrive at a lower bound

[TABLE]

For models with topologically deformable $U_{i}$ operators, like the 2D, 3D toric code models, our lower bound is invariant under perturbation. Because we can always move the excitations deep inside $D$, such that the distance from any excitation to the boundary $d_{e}\gg\xi\sim R$. For large subsystems, we would have $M=M^{\prime}$ and we do not lose any $U_{i}$, $U_{i}^{def}$ and $W_{i}$.

For models with $U_{i}$ not topologically deformable, e.g. the X-cube model and fractal spin liquids. There is usually some excitation that could not be moved deep inside subsystem $D$. Therefore, after adding perturbations, we typically lose a few $U_{i}$ and $U_{i}^{def}$ , such that $\bar{M}<M$ and $\bar{N}<N$. But for large subsystems which possess extensive $S_{nonlocal}$ before perturbation is added, this modification is small comparing to the leading contribution.

To summarize, for local perturbations satisfying assumption $\mathbf{S}$, and subsystem sizes much larger than $\xi$ and $\xi_{\alpha}$ we expect:

[TABLE]

For conventional topological orders $\mu=0$, and for fracton topological orders $\mu>0$ being a number depends on the model and subsystem geometry.

V Discussion and Outlook

In this paper, we have obtained a lower bound of the nonlocal entanglement entropy $S_{nonlocal}$ from assumptions about the topological excitations and the ground state condensates of Abelian topological orders and applied our method to several examples. For conventional topological orders, e.g. the 2D toric code model and the 3D toric code model, our lower bounds are saturated and topologically invariant. Whenever the lower bound is saturated, we get an explicit construction of a conditionally independent density matrix $\sigma^{\ast}$. For fracton topological orders Chamon (2005); Bravyi et al. (2011); Haah (2011); Yoshida (2013); Vijay et al. (2015, 2016), e.g. the X-cube model and the Sierpinski prism model, our lower bound depends on the geometry of the subsystems and $S_{nonlocal}$ is extensive for certain subsystem choices.

This method observes an intimate relation between $S_{nonlocal}$ and the topological excitations and the ground state condensates, and it obtains a lower bound of $S_{nonlocal}$ without calculating the entanglement entropy of any subsystem. A nonzero lower bound of $S_{nonlocal}$ is a result of the nonlocal nature of topological excitations, i.e. the fact that topological excitations could not be created alone by local operators. This nonlocal nature of topological excitations does not guarantee the operators which create the topological excitations to be topologically deformable and $S_{nonlocal}$ is not necessarily a topological invariance. Geometry-dependent $S_{nonlocal}$ is what appears in fracton models. It is beyond an established paradigm, i.e., the topological entanglement entropy, and should be treated as its generalization. The stability of the lower bound is discussed for local perturbations satisfying assumption $\mathbf{S}$, which should cover a large class of interesting local perturbations.

The different behaviors of $S_{nonlocal}$ may be used to distinguish fracton topological orders from conventional topological orders. Together with other methods being developed so far Kim (2013); Haah (2014), our result provides a better understanding of the entanglement properties of fractal models. Furthermore, the lower bound suggests (but not prove) different behaviors of $S_{nonlocal}$ between type I and type II fracton models. These different behaviors may be proven or disproven by later works.

Some of the assumptions in our method do not apply to non-Abelian models, a variant of our lower bound of $S_{nonlocal}$ for non-Abelian models is presented in Shi and Lu (2018). Also, it might be interesting to investigate possible implications of our method on relations among topological order, topological entanglement entropy and quantum black holes McGough and Verlinde (2013); Rasmussen and Jermyn (2017).

Acknowledgement

B.S. would like to thank Fuyan Lu for a comment on one assumption, Jeongwan Haah for providing a reference about Bravyi’s counterexample, Isaac H. Kim for a discussion and sharing one conjecture, and Michael Levin for a discussion about the entanglement in non-Abelian phases. This work is supported by the startup funds at OSU and the National Science Foundation under Grant No. NSF DMR-1653769 (BS,YML).