Search

Collections

Journals

About

Contact us

My IOPscience

Entanglement entropy and entanglement spectrum of triplet topological superconductors

This content has been downloaded from IOPscience. Please scroll down to see the full text. 2014 J. Phys.: Condens. Matter 26 425702 (http://iopscience.iop.org/0953-8984/26/42/425702) View the table of contents for this issue, or go to the journal homepage for more

Download details: IP Address: 152.3.102.242 This content was downloaded on 06/10/2014 at 11:25

Please note that terms and conditions apply.

Journal of Physics: Condensed Matter J. Phys.: Condens. Matter 26 (2014) 425702 (11pp)

doi:10.1088/0953-8984/26/42/425702

Entanglement entropy and entanglement spectrum of triplet topological superconductors T P Oliveira1,2 , P Ribeiro1,3,4 and P D Sacramento1,5 1

Centro de F´ısica das Interacc¸o˜ es Fundamentais, Instituto Superior T´ecnico, Universidade de Lisboa, Av. Rovisco Pais, 1049-001 Lisboa, Portugal 2 Departamento de F´ısica, Universidade Federal de Santa Catarina, 88040-900 Florian´opolis, Santa Catarina, Brazil 3 Max Planck Institute for the Physics of Complex Systems, N¨othnitzer Str. 38, 01187 Dresden, Germany 4 Max Planck Institute for Chemical Physics of Solids, N¨othnitzer Str. 40, D-01187 Dresden, Germany 5 Beijing Computational Science Research Center, Beijing 100084, People’s Republic of China E-mail: [email protected] Received 2 May 2014, revised 3 August 2014 Accepted for publication 20 August 2014 Published 2 October 2014 Abstract

We analyze the entanglement entropy properties of a 2D p-wave superconductor with Rashba spin–orbit coupling, which displays a rich phase-space that supports non-trivial topological phases, as the chemical potential and the Zeeman term are varied. We show that the entanglement entropy and its derivatives clearly signal the topological transitions and we ﬁnd numerical evidence that for this model the derivative with respect to the magnetization provides a sensible signature of each topological phase. Following the area law for the entanglement entropy, we systematically analyze the contributions that are proportional to or independent of the perimeter of the system, as a function of the Hamiltonian coupling constants and the geometry of the ﬁnite subsystem. For this model, we show that even though the topological entanglement entropy vanishes, it signals the topological phase transitions in a ﬁnite system. We also observe a relationship between a topological contribution to the entanglement entropy in a half-cylinder geometry and the number of edge states, and that the entanglement spectrum has robust modes associated with each edge state, as in other topological systems. Keywords: topological superconductors, entanglement entropy, topological entanglement entropy (Some ﬁgures may appear in colour only in the online journal)

1. Introduction

put forward [2]. For those, the characterization and detection of topological phases may be achieved by certain topological invariants associated with the ﬁlled energy bands. In the case of certain insulator classes, such invariants can be associated with direct physical response functions such as the quantized Hall conductivity [3]. On the experimental side, an increasing activity has been seen in the search for new topological insulators and topological superconductors. Various experimental signatures have been proposed and considerable experimental evidence has by now been found [4–21].

Topological phases of matter are characterized by global entanglement and correlations. Due to their non-local nature, topologically ordered phases are robust to local perturbations and have therefore received wide attention in the context of error-free quantum computation and quantum information processing [1]. Although a full classiﬁcation of topological phases is far from being achieved, very important steps have been given in this direction, in particular in the context of gapped noninteracting systems where a complete classiﬁcation has been 0953-8984/14/425702+11$33.00

1

© 2014 IOP Publishing Ltd

Printed in the UK

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

The absence of a local order parameter and the impossibility to use Landau symmetry breaking arguments turn topological phases into one of the most successful examples where entanglement measures, such as the entanglement entropy, bring new insights that could be hardly achieved by other more traditional methods. The use of quantum information concepts in condensed matter physics has already been extensively considered with entanglement measures being used to probe properties of many body states [22] and to detect quantum phase transitions [23]. Kosterlitz–Thouless transitions, lacking a local order parameter, were successfully detected calculating the ﬁdelity susceptibility of the XXZ spin chain [24, 25]. Many other systems were also studied such as the 1D Hubbard model [26, 27], spin-1/2 particles on a torus [28], the toric code model and the quantum eight-vertex model [29], the spin honeycomb Kitaev model [30–34], as well as other spin systems [35, 36]. In the thermodynamic limit, the ground state of a local Hamiltonian with a ﬁnite energy gap is characterized by shortrange correlations of local observables. For such systems, the ground state entanglement entropy of a subsystem A: SA = TrρA ln ρA (with ρA the reduced density matrix of the subsystem) follows the so-called area law [37], i.e. it is proportional to the size of the boundary of A in the limit of asymptotically large systems. Nonetheless, sub-leading corrections are expected. For a 2D system with perimeter P , the area law [37] translates to: SA = ξA P − γA + · · · .

contributions that may add up to the topological ones [45]. These corrections are negative and are given by the logarithm of the number of degenerate groundstates [45]. For systems with gapped quasiparticle-like excitations, but with gapless collective modes resulting from a spontaneous breaking of a continuous symmetry, there are also negative corrections that diverge logarithmically as the system size grows [46]. It is known that in certain systems the TEE vanishes [42], even for systems that cannot be adiabatically connected to trivial band insulators and there is no topological order in the sense of [47]. However, there are nontrivial topological phases as indicated by ﬁnite Chern numbers or winding numbers and edge states, such as in topological insulators or topological superconductors. These non-interacting fermionic systems have symmetry-protected topological properties (as evidenced by the protected gapless edge states or degenerate boundary states) but only have short-range entanglement [48, 49] in contrast to the long-range entanglement characteristic of the systems that have topological order. A way to differentiate the two types of entanglement is by noting that the shortrange entanglement systems may be reduced to trivial systems by the application of local unitary transformations [50], in contrast to the long-range entangled systems with topological order. A way to obtain topological order out of a noninteracting system is by the introduction of some projection operator, such as a Gutzwiller like transformation, as in the mapping of a p + ip spinless superconductor to the Moore– Read topologically ordered state [51]. For a generic domain A, γA may be non-zero and depends on the particular domain geometry. This has been observed for instance in [41] for a px + ipy spinless superconductor for which γA is ﬁnite and proportional to the number of corners of the partition. Even though the system has both a trivial and a nontrivial topological phase with quantum dimension D = 2, γTopo vanishes in both and is not able to distinguish between these phases. A similar situation was identiﬁed in the context of (the quadratic) Kitaev’s model where for a cylinder geometry an extra contribution to the entanglement entropy √ (ln 2) was identiﬁed [52] indicating that in some cases γTopo is not appropriate to identify and characterize topological phases. In this work, we calculate the entanglement entropy of a p-wave superconductor where the presence of Rashba spin– orbit coupling and a Zeeman term allow for a phase diagram that spans several topological nontrivial phases. We show that the entanglement entropy and appropriate derivatives of the entanglement entropy signal the topological transitions, as in other quantum phase transitions, topological or nontopological. We ﬁnd numerical evidence that the derivative of the entanglement entropy with respect to the magnetization seems to be related to the Chern numbers of each phase and displays approximate plateaus that scale in a way related to the relative Chern numbers. However, this approximate behavior does not seem to be universal and may be model dependent. We analyze both the extensive and non-extensive contributions in different geometries: a subsystem of an inﬁnite system and half of a cylinder. We also show that even though the TEE (γTopo ) vanishes everywhere in the thermodynamic limit, a suitable choice of geometry enables the identiﬁcation of the various

(1)

Here ξA is a non-universal constant term, γA contains the nonextensive contributions and · · · denote other contributions that vanish as P → ∞. If a phase has topological order, where some long-range entanglement is expected, then there should be an entropy reduction with respect to the area law of a nontopological gapped system. The non-extensive corrections have been previously shown to encode subtle effects due to the presence or absence of topological and/or long-range order [38]. Some contributions to γA are universal, in the sense that they are robust to transformations that do not close the gap and can be assigned to a given phase of mater [39, 40]. However, in some cases this information can be masked by non-universal terms [41, 42]. A prescription to extract the topological content of the entanglement entropy was given by Kitaev and Preskil [39] and Levin and Wen [40]. In their proposals they engineer particular ways to extract a universal topological contribution, here denoted as γTopo , from the non-extensive terms. This quantity, dubbed topological entanglement entropy (TEE), has been proposed to characterize certain topological phases and is related to the quantum dimension D of the system: γTopo = ln(D) [39, 40]. Since then, γTopo has been used as a signature of topological order in several systems as, for instance, in frustrated quantum dimer models and in the Kitaev honeycomb model [28, 43], as well as applied to detect topological order in spin liquid states [38, 44]. A rich set of different behaviours have been found, in particular for gapped systems that spontaneously break a discrete symmetry, where the correction to the area law has been observed to receive other 2

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

topological phases, in agreement with other signatures such as the Chern and winding numbers in a way similar to other topological systems. In section 2, we present the model studied and discuss its properties. In section 3, we present the formalism to calculate the entanglement entropy and the entanglement spectrum of a quadratic system. In section 4, we present the results. In section 4.1, we present the results for the entanglement entropy and its derivatives and their relation to the phase diagram determined using the Chern numbers of each phase. In section 4.1.1, we analyze the extensive contribution to the entanglement entropy and recover the overall characteristics of the phase diagram calculating the magnetization derivative of the entanglement entropy. In section 4.1.2, we consider the non-extensive contributions of a subregion of an inﬁnite system and compare the results as a function of the shape of the various subregions. In section 4.2, we show that, even if γTopo vanishes for an asymptotically large subregion, ﬁnite size effects are considerable and can be seen to clearly signal phase transitions. In section 4.3, we analyze the case of a cylinder geometry and show that the non-extensive contribution to the entanglement entropy can be related to the number of edge states. We conclude with section 5.

Figure 1. Topological phases and their Chern numbers, C, as a function of chemical potential, εF , and magnetization, Mz . γC is the TEE on a cylinder geometry (see text). The phases with C = 0 and γC = 0 are topologically trivial gapped phases.

The Pauli matrices σx , σy , σz act on the spin sector, and σ0 is the identity. We consider a unitary pairing contribution that reads Δ = i (d · σ ) σy with the vector d = (dx , dy , dz ) specifying the particular p-wave superconducting pairing. As the spin– orbit term breaks parity, a singlet pairing contribution is in principle also allowed but will not be considered here, for simplicity. Furthermore, we concentrate on the strong spin– orbit limit where the spin–orbit coupling is expected to be aligned with the pairing vector d.s = ||s|| ||d|| [57]. Even though in each topological phase the spin–orbit coupling may be turned off without affecting the topology, its presence, together with the TRS breaking Zeeman term, leads to a nonvanishing anomalous Hall effect and a ﬁnite Hall conductivity that may be used to obtain information about the topological phases [21]. In the absence of the Zeeman ﬁeld (Mz = 0), TRS is preserved and the system belongs to the symmetry class DIII. For this class the topological invariant belongs to Z2 [2, 58, 59]. In the presence of a Zeeman term, TRS is broken and the system belongs to the symmetry class D. The topological invariant that characterizes this class phase is the ﬁrst Chern number C, and the system is said to be a Z topological superconductor. Figure 1 shows the phase diagram, labeled by the Chern number values of each phase, as a function of chemical potential, εF , and Zeeman ﬁeld, Mz .

2. Triplet topological superconductor

Superconductivity with non-trivial topology may be obtained in several different ways [5]. However, non-centrosymmetric superconductors such as CePt3 Si [53]or the more recently discovered Li2 Pdx Pt3−x B [54], where s and p-wave paring may coexist, constitute rather natural candidates from an experimental viewpoint. In this work, we consider a model recently proposed in [55] that includes a Rashba-type spin– orbit coupling leading to an admixture of s- and p-wave pairing and Zeeman splitting terms that breaks time-reversal symmetry (TRS). For simplicity we consider only a triplet paring term with p-wave symmetry. In the presence of a TRS breaking Zeeman term and Rashba spin–orbit coupling, the phase diagram displays various transitions between topological phases characterized by different Chern numbers, as studied before in [21, 55, 56]. The Hamiltonian is given by 1 † Hˆ = C H (κ)Cκ , 2 κ κ where

H (κ) =

H0 (κ) Δ† (κ)

Δ(κ) −H0T (−κ)

3. Entanglement entropy and entanglement Hamiltonian

(2)

For quadratic Hamiltonians, the entanglement entropy can be obtained calculating the eigenvalues of the single particle correlation matrix deﬁned entirely in the subregion A [60–62]. In the following, we brieﬂy recall this procedure for a generic superconducting system and introduce some notation used in the subsequent sections. We use A¯ to denote the complement of A. A quadratic many-body Hamiltonian can be written in the form H = 1/2C † HC where H = H † is the one T † particle Hamiltonian and C = c1 , . . . , cN , c1† , . . . , cN is a column vector of annihilation and creation operators of N fermionic modes. A thermal density matrix of the composite

,

(3)

T † † , c and κ is the wave vector in the Cκ = cκ↑ , cκ↓ , c− κ↑ −κ↓ xy plane (the lattice constant is set to unity). The normal state Hamiltonian is given by H0 (κ) = κ σ0 − Mz σz + HR (κ) where κ = −2t (cos kx + cos ky ) − F is the kinetic part, with t the hopping amplitude and F the chemical potential; Mz is the Zeeman splitting ﬁeld (in units of energy); and H0 (κ) = s · σ is the Rashba spin–orbit term with s = α(sin ky , − sin kx , 0). 3

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

¯ is thus given by ρ = e− 2 C Ω C /Z, with system A + A, 1 † Ω = β H , β the inverse temperature and Z = Tre− 2 C Ω C =

1/2 . For a electronic system whose density det 1 + e−Ω matrix is of the form just described, the reduced density matrix of any subsystem A is itself of the quadratic form † ρA = e−1/2CA ΩA CA /ZA (with CA a vector of annihilation and creation operators restricted to subsystem A). The matrix ΩA can be most simply obtained from the correlations matrix † χi,j = Ci Cj as [60]: 1

†

ΩA = − ln χ−1 A −1

Figure 2. Domain A and its complement A¯ for the two considered geometries; left: subregion of an inﬁnite system; right: half of a section of a hollow ﬁnite cylinder.

where ν (ω) = α δ (ω − εα ) is the particle-hole symmetric [ν (ω) = ν (−ω)] density of states of single-particle entanglement Hamiltonian ΩA and −1 nf (ω) = 1 + eβω (9)

(4)

and ZA = det χ−1 A . χA;i,j is obtained from the global correlation matrix χ by restricting the indices i, j to label degrees of freedom of subsystem A only. As a result of this simpliﬁcation arising only for quadratic models, the entanglement entropy of a subsystem is given by 1 S [ρA ] = − (1 − λα ) ln (1 − λα ) + λα ln λα , (5) 2 α

is the Fermi function taken in equation (8) to have unit temperature β = 1. In this form the interpretation of equation (8) is particularly transparent: it is simply given by the sum of the entropies of the individual single particle states weighted by Fermi-statistics. Next, we study two different system geometries sketched in ﬁgure 2.

where λα are the eigenvalues of χA . The entanglement spectrum of a subsystem introduced by Li and Haldane [63] is deﬁned as the set of eigenvalues of the logarithm of the reduced density matrix: ΩA = − ln ρA + c, up to an overall additive normalization constant c. The many-body operator ΩA can be interpreted as an effective Hamiltonian if the subsystem A was taken to be in a Gibbs state at temperature T = 1. The entanglement spectrum was shown to contain information about states along the boundary between the two subsystems, including excited states [64–67]. When the topological phases result from inversion symmetry [68], it was shown to provide a more robust signature of the topology than the edge states. Its usefulness is also evident in the analysis of the spectral ﬂow of the entanglement spectrum, in particular the trace index, which may be related to changes in the Chern number in topological insulators [69, 70]. For quadratic systems, it is natural to deﬁne the so called ‘entanglement Hamiltonian’ [63] as ΩA = − ln ρA − ln ZA =

1 † C ΩA C . 2

3.1. Subregion of an infinite system

The correlation matrix of an inﬁnite translational invariant superconducting system can be written as d 2 k i k . (r − r ) χr,r = e χ (k ) , (10) (2π )2

with χ (k) = C k C †k = 1 − nf [H (k)], computed for a thermal Gibbs ensemble with temperature T = 1/β. For the ground-state one can simply take the β → ∞ limit, in which case nf becomes a Heaviside-theta function. Diagonalizing H (k) explicitly, we are thus able to compute χA by numerically performing the integral in equation (10) for r, r ∈ A. 3.2. Subregion of a finite cylindrical system

(6)

For a system placed in a cylindrical geometry (see ﬁgure 2), the correlations matrix χ(k) factorizes as a function of the momentum of the compact direction (here taken to be kx ):

The constant − ln ZA was chosen in order for the spectrum ofΩA to be particle-hole symmetric. The eigenvalues of ΩA can be written as Λn = (nα − 1/2)εα , (7)

χ r ,r =

1 ikx (x−x ) e χy,y (kx ) , Lx k

(11)

x

α:εα >0

where χ (kx ) = 1−nf [H (kx )] with H (kx ) the kx component of the Fourier-transformed Hamiltonian along the x direction, and with open boundary conditions along the y direction. Further restricting this matrix to a subsystem A that respects translational invariance along the x direction, one may deﬁne χA (kx ) = |y χy,y (kx ) y | (12)

where εα s are the eigenvalues of ΩA and the sum runs over non-negative εα . The vector n labels the eigenvectors of ΩA with its nα = 0, 1 being the occupation number of mode α. As ΩA is the single particle operator corresponding to ΩA , it is usually referred to as a single particle Hamiltonian and the set of εα s as the single particle entanglement spectrum [68–70]. Noting that εα = − ln(λ−1 α − 1) equation (5) can also be written S [ρA ] = − dω ν (ω) nf (ω) ln nf (ω) , (8)

y,y ∈A

as the kx -resolved single particle correlation matrix in the A domain. In the same way, the entanglement Hamiltonian 4

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

Figure 3. Left panel: local entanglement entropy of a squared L × L system with L = 6 as a function of chemical potential and

magnetization. Middle and right panels: derivatives of the entanglement entropy with respect to the chemical potential and the Zeeman term, respectively. The orange lines correspond to speciﬁc cuts in the phase diagram that are considered in the following.

can be kx -resolved: ΩA (kx ) = − ln χA (kx )−1 − 1 . As a function of kx , the entanglement entropy therefore factorizes S [ρA ] = kx Skx with Skx = I [εα (kx )] (13)

Even though in the superconducting phase a quantization is not expected, the plateau-like structure provides a sensible signature of the phases. As shown before [21], a similar result was obtained for the Hall conductivity and its derivatives. Even though the Hall conductivity is also not quantized, it allows a clear signature of the transitions between topological phases and, to some extent, provides information about the change in the Chern number, although it does not provide information about the actual Chern number. In that respect, the entanglement entropy seems to provide a more detailed information about the topological phases. This interesting approximate relation, between the derivative of the entanglement entropy and the Chern number does not seem to be universal, however and is model dependent. We have considered a model for a two-band 2D topological insulator, which displays phases with different Chern numbers, and we have calculated the entanglement entropy and its derivatives. The two-band model is of the type σ , where σ is the vector of Pauli matrices acting on the H = h· pseudospin space of the two bands, and the vector h has components hx = 2(cos kx +cos ky ), hy = 2(cos kx −cos ky ), hz = 2(sin kx + sin ky ) + 2(cos(2kx ) + cos(2ky )) + δ. Here, the parameter δ plays the role of a Dirac mass. As this parameter changes, the system changes its Chern number; as (δ < 0; 0 < δ < 4; 4 < δ < 8; δ > 8) the Chern number takes the values (C = 0; C = −1; C = 1; C = 0), respectively. In the regime of large δ, where C = 0, we ﬁnd that the derivative of the entanglement entropy with respect to δ is very small, as for the topological superconductor. The derivative signals the various transitions, as is observed in general quantum phase transitions. The derivative is ﬁnite and follows, approximately, a plateau like-structure when the system enters the phase with C = 1. However, as the phase with C = −1 is reached and a new plateau regime is reached, it turns out that, unlike in the topological superconductor, the derivative does not change sign as the Chern number changes from C = 1 to C = −1, as the results for the model considered in this paper might suggest. In the following, we study separately the extensive (proportional to the perimeter) and non-extensive contributions to the entanglement entropy, their dependence on the topological phase and on the geometry of the boundary. We ﬁrst consider the geometries displayed in the top panel of ﬁgure 5, then we compute the entanglement entropy using the method explained in section 3.1 and analyse the area law

α

the contribution of each momentum sector, where εα (kx ) are the eigenvalues of ΩA (kx ) and 1 1 1 1 1 I (ω) = − (14) ln + ln 2 1 + e−ω 1 + e−ω 1 + eω 1 + eω 4. Results 4.1. Entanglement entropy, derivatives and Chern numbers

We consider ﬁrst a subregion of an inﬁnite system. The entanglement entropy, S[ρA ], computed for a square L × L subregion of size L = 6, as a function of chemical potential and Zeeman term is given in ﬁgure 3, left panel, to illustrate the general behavior. The transition lines of ﬁgure 1, although visible, are not particularly well-deﬁned for such subsystem size. However, taking derivatives of the entanglement entropy with respect to the chemical potential or the Zeeman term, as shown in ﬁgure 3 (center and right panels), reveals that, even for such small sub-systems, the entanglement entropy clearly signals the transition lines displayed in ﬁgure 1. Note that the derivative with respect to the Zeeman term preserves the symmetry of the phase diagram as a function of the chemical potential. Figure 4 shows the convergence of the rescaled derivatives of the entanglement entropy with system size, along the three phase-space cuts of ﬁgure 3. We see that the convergence is fast even though we consider relatively small values of L. Away from the phase transition lines, where the entanglement entropy clearly signals the transition, we observe that the derivative of the entanglement entropy is mildly varying within each phase displaying a plateau-like structure. Indeed, inside each phase the derivatives reach a set of values that are approximately proportional to the value of the Chern number of the respective phase, with a proper rescaling of about 0.1. Particularly, in the regime where the Chern number vanishes the entanglement entropy is fairly independent of the magnetization. This has also been observed for other phase-space cuts. 5

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

Figure 4. Magnetization derivative of the local entanglement entropy of a square geometry with L × L sites for L = 3, 6, 10 as a function

of Zeeman ﬁeld Mz for F = −1, −3 and as a function of the chemical potential F for Mz = 0.5.

Figure 5. Top panel: different geometries of subsystems used for computing the entanglement entropy. The red dots correspond to the considered lattice degrees of freedom. The black lines deﬁne the boundary of a given geometry and were used to compute its perimeter P . Bottom panel: coefﬁcient of the perimeter contribution to the entanglement entropy, ξA , for different geometries, A = S, L, ..., obtained ﬁtting the numerically computed entanglement entropy for several system sizes and plotted as a function of Mz with F = −1. Inset: ratios of ξ ’s as a function of Mz .

for these systems. In the numerical calculations, we typically consider systems of linear size L 30.

rotated square (rS), rotated L-shaped (rL) and triangular-like shape (T), see the top panel of ﬁgure 5. We analyse their inﬂuence on the extensive contribution ξA to the entanglement entropy according to equation (1), where A assume the shapes S, L, U, rS, rL or T. Figure 5 (bottom panel) shows the various extensive contributions (ξS , ξL , ξU , ξrS , ξrL and ξT ) obtained by ﬁtting the

Entanglement entropy: extensive contribution. As discussed above, the geometry of the subregion directly affects the results of the entanglement entropy. Here we consider the following geometries: square (S), L-shaped (L), U-shaped (U),

4.1.1.

6

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

Figure 6. Left panel: sub-leading correction for the geometries L and U as a function of the correction of the S geometry. Middle panel: sub-leading correction for the geometries rL as a function of the rS one. Right panel: contribution of the different kinds of corners to the entanglement entropy.

asymptotic large perimeter (P ) limit for the different shapes along a cut in phase space where εF = −1 is kept constant. The general trend is a decrease of these contributions since as the magnetization increases, the entropy decreases. The S, L and U shapes yield the same contribution. The values for the rS and rL shapes also coincide but present a larger value than the one for S, L and U. This phenomenon occurs as the two kinds of boundaries, labeled i and ii in the ﬁgure, have different orientations with respect to the underlying lattice. The√T shape has mixed boundaries, having √ a fraction √ 2/(2 + 2) of type i boundary and a fraction 2/(2 + 2) of type ii. The values of ξT lay between ξS and ξrS . The inset shows that the ratios ξrS /ξS and ξT /ξS vary with Mz and are therefore model speciﬁc. Also, it is shown that the ˜ behavior √ of ξT with√Mz can be reproduced deﬁning ξT = (2ξT + 2ξrS )/(2 + 2). From these observations, we conclude that the extensive contribution varies with the speciﬁc details of the boundary and the ratio between different kinds of boundaries is not universal, depending on the speciﬁc details of the system. Nonetheless, when mixed boundaries are considered, their effect is additive and the total extensive contribution is an average over the ξ ’s of the different types of boundaries. 4.1.2.

γII when considering the type II corners. The numerical veriﬁcation of this fact enables us to compute the relative contribution of each kind of corner to the entanglement entropy. Figure 6 (right panel) shows the comparison of the different contributions as a function of Mz . We can see that there is a general decrease of these contributions as the magnetization increases, as expected, except at the transition points. The ratios γII /γI and γIII /γI are dependent on the details of the system and are therefore non-universal quantities. A similar conclusion arises if plotting γA as a function of ξA , for any subsystem A. Therefore, the results for the nonextensive contribution to the entanglement entropy do not show topological nature. Note that the relations γL = 3/2γS and γU = 2γS obtained previously imply the absence of a topological term, since if present this would add a constant term to these relations. 4.2. Vanishing of TEE

For a generic system, the term γA from equation (1) may have contributions coming from the corners of the geometry as well as a topological one. However, as mentioned above, this term for the model and geometries considered has no topological contribution. As discussed in the Introduction, this is expected for certain systems [41, 42] including quadratic Hamiltonians even with non-zero Chern numbers or Z2 phases due to the presence of short-range entanglement instead of the long-range entanglement present in systems with topological order. The absence of this long-range entanglement precludes the entropy reduction that is measured by the topological entanglement entropy. In order to conﬁrm the absence of the topological contribution, we consider the method described in [39], which allows us to compute the topological term. The procedure is designed to minimize ﬁnite size effects and isolate a possible topological contribution. Here, the so-called TEE is obtained considering three regions A, B, C (deﬁned in ﬁgure 3 of [39]) and calculating

Entanglement entropy: non-extensive contributions.

We now turn to the study of the non-extensive contributions (γS , γL , γU , γrS , γrL and γT ) to the entanglement entropy. γA for each geometry is obtained by ﬁtting the sub-leading term of the entanglement entropy. Figure 6 (left panel) shows the results for γL (γU ) and γS corresponding to the range of values [0,6] of the Zeeman coupling Mz at ﬁxed εF = −1. This plot clearly shows that all the non-extensive contribution is due to the corners of the geometry, since it does not distinguish between different topological phases, as observed before in other problems [41, 52, 71]. The effect of type I corners (see right panel) on the non-extensive contribution γI to the entanglement entropy is additive [41] as one obtains γS /4 = γL /6 = γU /8, i.e. γA = nA γI with nA the number of corners and A assuming the shapes S, L or U. A similar conclusion can be deducted from ﬁgure 6 (middle panel) to the non-extensive contribution

STopo = SA + SB + SC − SAB − SBC − SAC + SABC = − 2γTopo . 7

(15)

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

Figure 7. Left panel: γTopo for L = 24 as a function of chemical potential and magnetization. Middle and right panels: ﬁnite-size effects of

γTopo along cuts in the phase diagram.

In the ﬁrst panel of ﬁgure 7, we present results for L = 24, as a function of the chemical potential and Zeeman term. We consider the three regions A, B, C immersed in an inﬁnite system, for which we may calculate the correlations matrix in terms of the momentum space solution of the topological superconductor. Even though in the thermodynamic limit the TEE vanishes, for a ﬁnite system the results clearly shows the transition lines between the various phases. A close look into the numerical results shows that inside the various regions, the TEE is small but ﬁnite for ﬁnite systems. However, as shown in the other panels of ﬁgure 7, far from the transition lines, γtopo → 0, for all phases without distinction, as the system size grows. This is to be expected for a Z superconductor (related to the integer quantum Hall effect for which there is no topological order, in the sense of Wen [47, 72]). Therefore, in these geometries the TEE does not distinguish the various phases. However, in the next section, we will see that a ﬁnite system with cylindrical geometry allows us to isolate the term that has information about the topological nature of the various phases [38, 52]. This geometry has a different topological structure of the partitions considered in sections 4.1 and 4.2. The cylinder symmetry is equivalent to a disk shape where the inner region has a hole (as a coin with a hole). The subregion A and its complementary subregion are then partitions of the remaining circular system.

Concerning the contribution that is extensive in the perimeter of the system P = Lx , we ﬁnd, as expected, a very good agreement between ξS (see section 4.1.2) and ξCy , as in both cases the boundary is of type i (see ﬁgure 5). This geometry has no corners and the non-extensive component of the entanglement entropy, γCy , receives contribution from edges states only. We ﬁnd that in the topologically trivial phases [21, 55] with C = 0 (no edge states), γCy vanishes. In the non-trivial topological phases, we ﬁnd that γCy ∼ − ln 2 when C = 0 or C = −2 (both phases have four edge states) and γCy ∼ −1/2 ln 2 when C = 1 or C = −1 (with two edge states), as indicated in ﬁgure 1. Thus, in this case the contributions to the TEE can be written as: γCy = −(1/2 ln 2)gedge , where gedge is the number of √ edge state pairs. Thus, we identify each contribution as ln 2. This fractional entropy contribution is characteristic of a half-fermion (Majorana). A similar result has been found for the Majorana mode of the n-channel Kondo model with n = 2, S = 1/2, due to the impurity contribution [73]. To understand the origin of these terms, we analyze each α contribution to the entanglement entropy, as given in equation (5). All contributions to the entropy are positive deﬁnite. Therefore, for any system size the overall entropy is also positive. However, as the system size along the height of the cylinder changes, the total entropy may extrapolate to a positive or a negative value (at zero system size). In a system with long-range entanglement, the extrapolation is negative and the TEE is positive. In our case, it extrapolates to a positive value. We ﬁnd that there is a set of robust terms of values +(1/2) ln(2), that are independent of the system size (therefore, they are not affected by the extrapolation to zero system size). These terms are related to the edge states (see below). On the other hand, extrapolating the remaining nonrobust terms to the zero size limit, gives a negative contribution to the entanglement entropy. As it turns out, summing the contributions of the robust modes with the extrapolated value of the remaining contributions gives an overall positive entropy (negative value to γCy ). We have checked that in the ﬁnite system with a square geometry (open boundary conditions

4.3. Entanglement entropy in cylinder geometry

We now address the entanglement entropy on a cylindrical geometry (Cy) (see ﬁgure 2). The subregion A is in this case half of the (hollow ﬁnite) cylinder. For the numerical calculations, we considered cylinders of lengths up to Ly = 70 and perimeters of varying sizes up to Lx = 400. We consider open boundary conditions along the y direction and periodic boundary conditions along the x direction. In the picture of the disk with the hole in the middle, this means that the cut into two subregions implies that both regions have two edges. In the conﬁgurations considered in sections 4.1 and 4.2, the inner region only has the outer edge. 8

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

Figure 8. Upper panel: single-particle entanglement spectrum as a function of kx computed for a representative point of the different topological phases of ﬁgure 1, corresponding respectively to {F , Mz } = {1., 2.}, {2., 1.}, {6., 1.}, {4., 2.}, {2., 4.} (from left to right) with Lx = 400 and Ly = 70. The blue arrows point to extra states arising at εα = 0 that yield the extra 1/2 ln 2 to the entropy. Lower panel: contribution to the entanglement entropy of each momentum sector. The red points correspond to an extra contribution of 1/2 ln 2.

trivial states: for C = 0, small deformations of the entanglement ‘bands’ do not lead to a gapless entanglement spectrum. For C = 0 two situations may arise: if there are no edge modes the entanglement spectrum is gapped (see ﬁgure 8 third panel); in the presence of edge modes, it is still possible to open up a gap by lifting the degeneracy arising at kp (see ﬁgure 8 second panel with kp = −π ). The degenerate states are associated with the gapless states along the two edges of subregion A. The Majorana states could then be established at each edge and the nonlocal entanglement between the two would give rise to the two states that are shown in the entanglement spectrum. The entanglement between these two degenerate states leads to a fermion that may be either present or absent and explains the degeneracy in the entanglement spectrum. Since this entanglement is non-local, it does not depend on the size or geometry. In [69], a related but somehow different observation has been made by Hughes et al for the single-particle entanglement spectrum of Chern and spin Hallinsulators. In their work, the authors report that, for a trivial insulator with edge modes, even if the entanglement spectrum is gapless the bands that cross zero are disconnected from the rest of the spectrum.

in both spatial directions) there are no robust modes, i.e. all contributions to the entanglement entropy are not robust to changes in Lx . As will be shown next, these robust contributions are associated with the edge states. Entanglement spectrum. In order to further understand these results, we analyze in detail the entanglement spectrum for the cylinder geometry. Figure 8 shows the single-particle entanglement spectrum for each of the phases of ﬁgure 1 as a function of kx . Note that the Hamiltonian is particle-hole symmetric Ω (kx ) = −Ω (−kx ). When the entanglement spectrum is gapless, the zero energy mode is doubly degenerate. This feature has been observed previously in a superconductor [65] and in the study of the entanglement spectrum of the Kitaev model [74]. As a consequence, in the calculation of Skx , the degenerate states are responsible for an extra 1/2 ln 2 contribution: Skx = s˜ (kx ) + 1/2 ln 2 p δkx ,kp , with s˜ (kx ) a continuous function of kx and where the kp ’s (= 0, −π in the examples of ﬁgure 8) are the values of the momentum for which the degeneracy arises. The function s˜ (kx ) is plotted in ﬁgure 8 where the kp points are also identiﬁed. For large values of Lx , the entanglement entropy can be approximated by dkx S [ρA ] s˜ (kx ) Lx + (1/2 ln 2) gedge , (16) 2π

4.3.1.

5. Discussion

In this work, we have analysed the entanglement entropy of a 2D topological p-wave superconductor. The main conclusions of our work are the following: The entanglement entropy clearly signals the topological transitions and its derivatives have sharp features around the transition lines, even for small system sizes. Moreover, the derivative of the entanglement entropy with respect to the Zeeman term shows approximate plateaus that provide a sensible signature of each topological phase and seems to be connected to the Chern number. However, the result may be model dependent. We have checked that similar results may be obtained for topological insulators, but we have not been able to prove with generality the relationship. We separately analyzed the contributions to the entanglement entropy that are proportional to or independent of the

replacing the sum over momentum values by an integral. This x corresponds to ξCy = dk and γCy = − (1/2 ln 2) gedge . s˜ 2π (kx ) The gapless points are characterized by εα = 0 and, consequently, λα = 1/2. This translates into a contribution to the entropy of 1/2 ln(2). The entropy for each transverse momentum is given by a sum over all bands shown in ﬁgure 8. Away from the gapless points, the entropy is given by s˜ (kx ). As one approaches the gapless points, one of the contributions approaches the value 1/2 ln(2). At the gapless point, an extra state appears that is degenerate with the state that approaches the value 1/2 ln(2) and there is a jump of the entropy at that momentum value, signalled by the red points in the lower panels of ﬁgure 8. Figure 8 shows a striking difference between the entanglement spectrum of topological insulators and topologically 9

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

perimeter of the subsystem, as a function of the Hamiltonian coupling constants and the geometry of the ﬁnite subsystem. We have conﬁrmed that the extensive contribution, ξA , is generically non-universal depending on the speciﬁc details of the system. On the other hand, the non-extensive contribution, γA , is directly related to the corners of the geometry. Its relation is due to the orientation of the corners with respect to the underlying lattice, and smaller angles give a larger contribution to the entanglement entropy. The ratios of the non-extensive contributions between the various geometries show that some universality may be obtained for other lattice geometries and models. Having a way to estimate ξA based on the characteristic lengths of the problem would be insightful; however, in this study we found no clear relation of ξA with the characteristic decay length of various correlation functions. Due to the quadratic nature of the BCS Hamiltonian, it is expected that the TEE vanishes inside each topological phase. The present work veriﬁed this expectation but also showed that the TEE for a ﬁnite system still signals the topological phase transitions. In a cylinder geometry, the topological contribution to the entanglement entropy is ﬁnite and is in general negative [38]. Even though the system has no topological order, the gapless edge states characteristic of the symmetry-protected topology leave traces in the non-extensive contribution to the entanglement entropy. The reason is associated with a robust negative contribution due to the gapless edge states, each contributing (1/2) ln 2 to the entropy, as is characteristic of Majorana (fractional) edge states. The analysis of the entanglement spectrum [52] shows that the extra contributions to the entropy are associated with spectral degeneracies, in agreement with previous works for other topological systems. As a ﬁnal note, we stress that even though the TEE vanishes in this problem (in the thermodynamic limit and inside each region in the phase diagram) and shows no topological order in the sense of Wen, complementary information is obtained choosing appropriate boundary conditions. Detailed study uncovers information contained in the entanglement entropy and the entanglement spectrum that can be directly related to the topological phases, alternatively characterized by the presence of non-vanishing Chern numbers or the presence of edge states of Majorana type.

[3] Thouless D J, Kohmoto M, Nightingale M P and den Nijs M 1982 Phys. Rev. Lett. 49 405 [4] Hasan M Z and Kane C L 2010 Rev. Mod. Phys. 82 3045 [5] Qi X-L and Zhang S-C 2011 Rev. Mod. Phys. 83 1057 [6] Alicea J 2012 Rep. Prog. Phys. 75 076501 [7] K¨onig M, Wiedmann S, Br¨une C, Roth A, Buhmann H, Molenkamp L W, Qi X L and Zhang S C 2007 Science 318 766 [8] Hsieh D, Qian D, Wray L, Xia Y, Hor Y S, Cava R J and Hasan M Z 2008 Nature 452 970 [9] Xia Y et al 2009 Nat. Phys. 5 398 [10] Mourik V, Zuo K, Frolov S M, Plissard S R, Bakkers E P A M and Kouwenhoven L P 2012 Science 336 1003 [11] Das A, Ronen Y, Most Y, Oreg Y, Heiblum M and Shtrikman H 2012 Nat. Phys. 8 887 [12] Deng M T, Yu C L, Huang G Y, Larsson M, Caroff P and Xu H Q 2012 Nano Lett. 12 6414 [13] Zhang F, Kane C L and Mele E J 2013 Phys. Rev. Lett. 111 056402 [14] Nadj-Perge S, Drozdov I K, Bernevig B A and Yazdani A 2013 Phys. Rev. B 88 020407 [15] Veldhorst M et al 2012 Nat. Mater. 11 417 [16] Williams J R, Bestwick A J, Gallagher P, Hong S S, Cui Y, Bleich A S, Analytis J G, Fisher I R and Goldhaber-Gordon D 2012 Phys. Rev. Lett. 109 056803 [17] Rokhinson L P, Liu X and Furdyna J K 2012 Nat. Phys. 8 795 [18] San-Jos´e P, Cayao J, Prada E and Aguado R 2013 New J. Phys. 15 075019 [19] Schnyder A P, Brydon P M R and Timm C 2012 Phys. Rev. B 85 024522 [20] Ojanen T and Kitagawa T 2013 Phys. Rev. B 87 014512 [21] Sacramento P D, Ara´ujo M A N and Castro E V 2014 Europhys. Lett. 105 37011 [22] Amico L, Fazio R, Osterloh A and Vedral V 2008 Rev. Mod. Phys. 80 517 [23] Gu S J 2010 Int. J. Mod. Phys. B 24 4371 [24] Yang M-F 2007 Phys. Rev. B 76 180403 [25] Chen S, Wang L, Hao Y and Wang Y 2008 Phys. Rev. A 77 032111 [26] Gu S-J, Kwok H-M, Ning W-Q and Lin H-Q 2008 Phys. Rev. B 77 245109 [27] Campus Venuti L, Cozzini M, Buonsante P, Massel F, Bray-Ali N and Zanardi P 2008 Phys. Rev. B 78 115410 [28] Hamma A, Zhang W, Haas S and Lidar D A 2008 Phys. Rev. B 77 155111 [29] Abasto D F, Hamma A and Zanardi P 2008 Phys. Rev. A 78 010301 [30] Yang S, Gu S-J, Sun C-P and Lin H-Q 2008 Phys. Rev. A 78 012304 [31] Abasto D F and Zanardi P 2009 Phys. Rev. A 79 012321 [32] Trebst S et al 2007 Phys. Rev. Lett. 98 070602 [33] Zhao J-H and Zhou H-Q 2009 Phys. Rev. B 80 014403 [34] Wang Z, Ma T, Gu S-J and Lin H-Q 2010 Phys. Rev. A 81 062350 [35] Castelnovo C and Chamon C 2008 Phys. Rev. B 77 054433 [36] Eriksson E and Johannesson H 2009 Phys. Rev. A 79 060301 [37] Eisert J, Cramer M and Plenio M B 2010 Rev. Mod. Phys. 82 277 [38] Jiang H C, Wang Z and Balents L 2012 Nat. Phys. 8 902 [39] Kitaev A and Preskill J 2006 Phys. Rev. Lett. 96 110404 [40] Levin M and Wen X G 2006 Phys. Rev. Lett. 96 110405 [41] Ali N B, Ding L and Haas S 2009 Phys. Rev. B 80 180504 [42] Nussinov Z and Ortiz G 2009 Ann. Phys. 324 977 [43] Furukawa S and Misguich G 2007 Phys. Rev. B 75 214407 [44] Depenbrock S, McCulloch I P and Schollw¨ock U 2012 Phys. Rev. Lett. 109 067201 [45] Jiang H C, Yao H and Balents L 2012 Phys. Rev. B 86 024424

Acknowledgments

We thank discussions with V´ıtor Rocha Vieira and the partial support of the Portuguese FCT under grant PESTOE/FIS/UI0091/2011. PR was supported by the Marie Curie International Reintegration Grant PIRG07-GA-2010-268172. TPO would like to thank the Brazilian agency Coordenac¸a˜ o de Aperfeic¸oamento de Pessoal de N´ıvel Superior (CAPES) for the ﬁnancial aid. References [1] Nayak C, Simon S H, Stern A, Freedman M and Das Sarma S 2008 Rev. Mod. Phys. 80 1083 [2] Ryu S, Schnyder A P, Furusaki A and Ludwig A W W 2009 New J. Phys. 12 065010 10

J. Phys.: Condens. Matter 26 (2014) 425702

T P Oliveira et al

[46] Kallin A B, Hastings M B, Melke R G and Singh R R P 2011 Phys. Rev. B 84 165134 [47] Wen X-G 1995 Adv. Phys. 44 405 [48] Chen X, Liu Z-X and Wen X-G 2011 Phys. Rev. B 84 235141 [49] Lu Y-M and Vishwanath A 2012 Phys. Rev. B 86 125119 [50] Chen X, Gu Z-C and Wen X-G 2010 Phys. Rev. B 82 155138 [51] Read N and Green D 2000 Phys. Rev. B 61 10267 [52] Yao H and Qi X-L 2010 Phys. Rev. Lett. 105 080501 [53] Bauer E, Hilscher G, Michor H, Paul C H, Scheidt E W, Gribanov A, Seropegin Y, Noel H, Sigrist M and Rogl P 2004 Phys. Rev. Lett. 92 027003 [54] Togano K, Badica P, Nakamori Y, Orimo S, Takeya H and Hirata K 2004 Phys. Rev. Lett. 93 247004 Badica P, Kondo T and Togano K 2005 J. Phys. Soc. Japan 74 1014 [55] Sato M and Fujimoto S 2009 Phys. Rev. B 79 094504 [56] Oliveira T P and Sacramento P D 2014 Phys. Rev. B 89 094512 [57] Frigeri P A, Agterberg D F, Koga A and Sigrist M 2004 Phys. Rev. Lett. 92 097001 [58] Schnyder A P, Ryu S, Furusaki A and Ludwig A W W 2008 Phys. Rev. B 78 195125 [59] Schnyder A P, Ryu S, Furusaki A and Ludwig A W W 2009 Advances in Theoretical Physics ed V Lebedev and M Feigel’man (AIP Conf. Proc.) vol 1134 (Melville, NY: American Institute of Physics) p 10

[60] Peschel I 2003 J. Phys. A: Math. Gen. 36 L205 [61] Chung M C and Peschel I 2001 Phys. Rev. B 64 064412 Cheong S A and Henley C L 2004 Phys. Rev. B 69 075111 [62] Peschel I and Eisler V 2009 J. Phys. A 42 504003 Peschel I 2012 Braz. J. Phys. 42 267 [63] Li H and Haldane F D M 2008 Phys. Rev. Lett. 101 010504 [64] Regnault N, Bernevig B A and Haldane F D M 2009 Phys. Rev. Lett. 103 016801 [65] Fidkowski L 2010 Phys. Rev. Lett. 104 130502 [66] Poilblanc D 2010 Phys. Rev. Lett. 105 077202 [67] Sterdyniak A, Regnault N and Bernevig B A 2011 Phys. Rev. Lett. 106 100405 [68] Turner A M, Zhang Y and Vishwanath A 2010 Phys. Rev. B 82 241102 [69] Hughes T L, Prodan E and Andrei Bernevig B 2011 Phys. Rev. B 83 245132 [70] Alexandradinata A, Hughes T L and Andrei Bernevig B 2011 Phys. Rev. B 84 195103 [71] Papanikolaou S, Raman K S and Fradkin E 2007 Phys. Rev. B 76 224421 [72] Zhang Y, Grover T and Vishwanath A 2011 Phys. Rev. B 84 075128 [73] Emery V J and Kivelson S 1992 Phys. Rev. B 46 10812 [74] Yao H and Qi X-L 2010 Phys. Rev. Lett. 105 080501

11