Today’s date: July 12, 2019
Entanglement entropy and the complex plane of replicas
Ferdinando Gliozzi and Luca Tagliacozzo
Dip. di Fisica Teorica and INFN, Università di Torino,
Via P. Giuria 1, 10125 Torino, Italy
School of Physical Sciences, the University of Queensland, QLD 4072, Australia
The entanglement entropy of a subsystem of a quantum system is expressed, in the replica method, through analytic continuation with respect to of the trace of the power of the reduced density matrix . We study the analytic properties of this quantity as a function of in some quantum critical Ising-like models in 1+1 and 2+1 dimensions. Although we find no true singularities for , there is a threshold value of close to 2 which separates two very different ‘phases’. The region with larger is characterised by rapidly convergent Taylor expansions and is very smooth. The region with smaller has a very rich and varied structure in the complex plane and is characterised by Taylor coefficients which instead of being monotone decreasing, have a maximum growing with the size of the subsystem. Finite truncations of the Taylor expansion in this region lead to increasingly poor approximations of . The computation of the entanglement entropy from the knowledge of for positive integer becomes extremely difficult particularly in spatial dimensions larger than one, where one cannot use conformal field theory as a guidance in the extrapolations to .
The entanglement entropy has become a privileged measure of bipartite entanglement for pure states, that is of the entanglement between a subregion of the system and the rest of it. In order to compute it, one divides the system in the state , into two complementary subsystems and . The reduced density matrix of the region is obtained by tracing out the degrees of freedom inside the region , . The correlations between and give rise to a mixed state for whose von Neumann entropy is the entanglement entropy of the region (see e.g.  for reviews). Other quantities such as the Rényi ()  and the Tsallis entropies () 
are strictly related to the entanglement entropy of the block . Indeed both of them coincide with it in the limit .
Studying directly the entanglement entropy in the ground state of a chosen Hamiltonian requires the ability to obtain such ground state. The complexity of this, in general, grows exponentially with the size of the system (recently new methods have been proposed to approximate the ground state of two dimensional systems that only scale polynomially with the system size [4, 5]).
However, when in Eq. (1.1) is a positive integer , the Tsallis and Rényi entropies can be computed in the context of quantum field theory (QFT), without the explicit knowledge of the ground state. The procedure is known as the “replica method” [6, 7] and is based on the fact that the expression can be written as a partition function of the model on a sheeted Riemann surface (or its generalization in dimensions), divided by , where is the partition function of the original system . The entanglement entropy is obtained from it through an analytical continuation from positive integer to real using . Even if all values of were provided for , analytic continuation to real might not be uniquely determined, as examples on disordered systems show . It is therefore desirable to have some information on the analytical properties of .
Here is an ultraviolet (UV) cutoff, is the central charge of the associated CFT and is a non-universal function of which is known only in some special cases. In some cases even the full analytic continuation has been performed  allowing one to extract also the non universal additive corrections to the scaling of the entanglement entropy .
In , however, CFT results are not available and one thus needs to carefully study the analytical properties of before trying to perform any analytic continuations of it.
The purpose of this work is to perform such a study by exploring the analytic properties of as a function of in the whole complex plane for a class of quantum Ising-like critical models in both 1+1 and 2+1 dimensions. By doing this we uncover a highly non-trivial analytical behaviour of for real in the interval . In all the models considered, we find a threshold value in the interval which divides into two ‘phases’. In the region with the quantity behaves very smoothly as a function of ; by expanding it in a Taylor series centred around any of the , we extract alternating series with monotone decreasing coefficients, then, any truncation can be used to approximate , the error of this approximation being less than . For this reason we call the region ‘perturbative’ region.
On the contrary, in the other ‘phase’, the Taylor expansions of with has coefficients whose moduli are no longer monotone decreasing. They show a characteristic peak whose height increases with the size of (see for instance Fig.1). In this regime any truncation of the Taylor series becomes an increasingly poor approximation of as the size of grows. Thus we name the region ‘non-perturbative’ region.
We shall see that the threshold value coincides with the point where the expression changes sign, that is
In the models we have considered it turns out that the two different ‘phases’ are separated by a smooth cross-over in both and . However we cannot exclude that, in higher dimensions, this cross-over turns into a real phase transition. In fact it has been recently observed in Ref.  that the Rényi entropies in the model undergo a phase transition as a function of for close to 3. This happens in the same range of in which the models we analyse develop the cross-over between the perturbative and the non-perturbative regions.
The above scenario has important consequences when trying to extract the entanglement entropy using the replica trick. On one side the absence of a real phase transition makes the analytic continuation (it involves) possible. On the other side, the different behaviour of the two ’phases’ implies that finding it can be really a non trivial task.
In fact, it is easy to find a simple function which accurately reproduces for integer . However, if we make the ‘obvious’ extension to real , by assuming that this functional form is also valid for real , we obtain an incorrect value of the entanglement entropy.
By exemplifying this we also have the opportunity to compare two different methods. One is based on a recent algorithm proposed in Ref.  which allows us to accurately evaluate the reduced density matrix of a two-dimensional quantum system defined on a lattice of small size (in this paper we generally refer to this method as the Hamiltonian approach). The other one is based on the replica method and applies a simple technique described in Ref.  to directly measure as the vacuum expectation value of a suitable observable in a Monte Carlo simulation, in any Euclidean lattice field theory.
The two methods agree on the Tsallis entropies for a wide range of integer , once the UV cutoff of the two approaches are matched through a suitable conversion factor (this by itself is a non trivial result that provides a strong confirmation of the validity of both approaches). However there is a mismatch in the results for the entanglement entropy. Its value directly computed in the Hamiltonian approach (the correct one) is different from the one obtained in the Monte Carlo approach from the naive continuation of the behaviour of the Tsallis entropies for integer to real . This is a consequence of the intricate landscape of close to . The scenario we have outlined is common to all critical quantum systems in one and two dimensions we have analysed.
The contents of this paper are organised in several sections. In section 2 we describe the quantum models we want to study with different numerical methods and, where available, analytical results. In section 3 we describe both the Monte Carlo and the variational techniques we use to analyse quantum systems in 2+1 dimensions. Section 4 and 5 present the main results of this paper. In Section 4 we outline the presence of two ‘phases’ for by considering the features of its Taylor expansions. In section 5 we look for the Lee-Yang zeros and other possible singularities in the complex plane of and we find no singularity related to the transition between these two regions. We unveil the presence of a smooth cross-over between the two regions. In Section 6 we discuss the effects of this scenario on entanglement entropy calculations. There we also present the comparison between the numerical results obtained with Monte Carlo and variational techniques. We conclude with a discussion of the results in Section 7.
2 The models
We analyzed the reduced density matrix of the ground state near and at the quantum critical point in a series of one-dimensional spin chain systems and their two-dimensional generalizations. In particular, we considered the XY Hamiltonian
where denote first neighbor sites, are Pauli matrices and , . The subsystem is a block of neighboring nodes.
We applied mainly the numerical methods described in [11, 13], but we used also analytical results. Analytical calculations exist for the critical XX chain , for the non-critical XY chain in a field  as well as for the Heisenberg chain .
Simplifications arise when mapping the quantum chains into classical two-dimensional classical spin systems . The quantum XY chain can be mapped, at least for a part of the parameter range, into a classical Ising model on a triangular lattice . In this way simple formulas can be easily obtained. For instance, if the subsystem is a block of spins of size much larger than the correlation length , embedded into an infinite chain with periodic boundary conditions at zero temperature in the ordered phase ((i.e. ) the reduced matrix is independent of and can be written in the simple form 
where the parameter which sets the scale is a known function of and [20, 21]. The factor of two in front of (2.2) arises from the asymptotic degeneracy of the ground state in the broken symmetry for . The exponent of two in (2.2) reflects instead the fact that the segment in a chain with periodic boundary conditions has two points of contact, i.e. two interfaces with the rest of the system.
we introduced the Dedekind function
because it has a simple behaviour under the modular transformation , namely
and this leads to rapidly converging products. Precisely we have
Since , the first exponential tells us, according to [6, 7], that the associated CFT in the critical limit is the free fermion model with central charge as expected, however the exact expression (2.7) allows to discuss the analytical structure of the reduced density matrix in the whole complex replica plane. This discussion is postponed to Section 5.
For the other critical models, to our knowledge there are no analytical results available for the reduced density matrices, and we have to rely on a set of different numerical techniques. For the generic one dimensional critical XY chains we use the expression of the reduced density matrices of an interval in terms of the free fermionic modes described in Ref.  (for more recent results see e. g. Ref. ) and numerically diagonalise it. In the case of finite chains at the critical point we use the one dimensional version of the variational technique described in Ref.  that employs a Tree Tensor Network as an Ansatz for the ground state of the system.
In two dimensions we analyse the Ising model in a transverse field and the XX model on finite tori with both the Monte Carlo technique described in Ref.  and the variational technique of Ref.  that uses again a Tree Tensor Network. In 2D the subsystem is one half of a torus bounded by two parallel straight lines.
3 Monte Carlo simulations and variational algorithms
In QFT the quantity can be written as the vacuum expectation value of a suitable observable defined on a larger system composed of decoupled copies of the original system.
The canonical partition function of a dimensional quantum system at inverse temperature can be computed in a standard way by doing the Euclidean functional integral in a -dimensional hyper-cubic lattice over fields periodic under . Therefore the system composed of independent replicas of the original system is described by the th power of :
where is a field configuration associated to the th copy and is the Euclidean action.
In the replica approach to entanglement the subsystem establishes a process of transferring information among the replicas through a specific coupling: the lattice links coming out of nodes of and directed into the (Euclidean) time direction cyclically connect the copy with the copy [6, 7]. Let us denote with the corresponding coupled action. It is easy to see that the quantity
has the desired property. In fact its vacuum expectation value in the system of independent copies of the original system is
is the partition function of the system in the coupled replicas.
This method can be simply implemented in Monte Carlo simulations. Here we apply it to the case in which the system is a periodic lattices of size , where the longer direction is the Euclidean time. We verified that the Euclidean time is large enough to ensure the absence of measurable finite temperature effects. The Monte Carlo simulations were performed with a particularly efficient implementation for the Ising and state Potts model described in  and with the standard action of the isotropic Ising model on a cubic lattice
All simulations were made at the critical coupling .
With this method we compute directly the Tsallis entropies for integer . A discussion of the results we obtained is postponed to Section 6.
Using the techniques of Ref.  we also extract an accurate variational approximation to the ground state with a Tree Tensor Network of the corresponding quantum model. This is the quantum Ising model on a square torus at its critical point obtained by setting the parameters of the Hamiltonian in Eq. 2.1 to =0 and . This approach is limited by the exponential growth of the Hilbert space with the system size. The benchmark calculations presented in Ref.  show that, nevertheless, it can be used safely at least for tori of dimension up to . It has the strong advantage that it allows us to extract the complete spectrum of the reduced density matrix for arbitrary blocks. It thus provides a way to directly compute for an arbitrary complex and also to directly compute the entanglement entropy. A detailed explanation on how to compute the reduced density matrices in this approach is contained in Ref. . We also postpone the discussion of the numerical results for both the Tsallis and entanglement entropy extracted with this method to Section 6.
4 The two regions of Tsallis entropies
As we discussed in the introduction the Taylor expansions of can be used to study the behaviour of the Tsallis and Rényi entropies as a function of and identify two very different regions. In order to exemplify this idea we have to express the Taylor expansion of about in a convenient form. Let be the non vanishing eigenvalues of . At first the normalisation of implies
Through the chain of identities
we define the Taylor coefficient . Since the eigenvalues of are smaller than 1 the Taylor expansion is an alternating series, i.e. . This is an important feature when discussing the possible location of the singularities in the complex plane or in estimating the error we make in truncations of the Taylor expansions, as we shall see later. From the spectrum of we can explicitly compute the above Taylor expansion.
As an example, the method described in Ref.  can be used to compute the reduced density matrix of a block of spins of an infinite critical Ising chain. From it, we compute the moduli of the first Taylor coefficients of and plot them in Fig.1 (a) for four different values of . It is interesting to notice that the plot shows a peak at some value of and that its height increases with .
A similar behaviour can be observed in a quantum two-dimensional critical Ising model on a torus as shown in Fig.1 (b). There we plot again the absolute value of the Taylor coefficients of the expansion of about for four different values of . The subsystem in this case is half a torus. The reduced density matrices are obtained with the techniques described in Ref. . By comparing the same expansion for different systems sizes, we appreciate how the peaks move with and, more importantly, their height increases with .
By varying the expansion centre of the Taylor expansions we can check the behaviour of the coefficients in other regions of the axis. The expression for the coefficients of the Taylor expansion about is
We repeat the same study presented in Fig. 1 in the region and observe that the peaks of Fig. 1 (a) and (b) start moving to left as we increase the value of to eventually disappear at a threshold value . In the region where the peaks are present we cannot obtain a good approximation of by truncating the power series (4.2) down to only few terms. Also, the convergence of the series (4.2) deteriorates with the size of the system, preventing the possibility to find a good truncation scheme independent of the system’s size. For this reason we call the region the non-perturbative region *** Starting from here, we will gradually drop the notation in favour of . The reader should not be confused. Indeed the Taylor coefficient for do not depend on but rather on , the centre of the expansion. However has the same support than and we can safely simplify the notation..
For the coefficients become monotone decreasing in modulus. As a consequence, in this region, the partial sum can be used to approximate with an error smaller than . In this case thus, we call the region the perturbative region.
The threshold value is defined as the minimal value of such that
that is, the Taylor coefficient are monotonically decreasing functions of . In all the models we have considered, as soon as the criterion in Eq. 4.4 is fulfilled by the second and third coefficients, that is, , it is also fulfilled by all the others pairs. Using this property we can define as the solution of Eq.(1.3) anticipated in the introduction. We can also use
the sign of the expression (1.3) defining , as a sort of order parameter distinguishing between the perturbative and non-perturbative ‘phases’. As expected, the value of the threshold depends on the size of the system, i.e. . By analysing this dependence we found that in our critical systems obeys a simple power law
where in the 1D critical Ising model and in the 2D critical Ising model. It is also interesting to notice that for all the models we have analysed (see Figs. 2 (a) and (b)).
5 The nature of the transition
Having uncovered the presence of two different regions for Tsallis (and Rényi) entropies as functions of , we now want to investigate whether these two regions are well separated phases divided by a transition point or if instead they are different regions of the same phase smoothly connected through a cross-over. Of course, these two scenarios have very different physical implications. In the case where the two regions are separated by a true phase transition for some real in the range , there would be no possibility to analytically continue the Tsallis entropy from one phase to the other. This, for example, would imply that the entanglement entropy could not be obtained via the replica method. This is unlikely to happen in one dimensional critical systems since the entanglement entropy has been already successfully computed with the replica trick in some 1+1 CFT . In two dimensions however, for the models we are considering, the entanglement entropy has not yet been obtained analytically and one can be concerned about if with the replica trick is a viable method to perform such a calculation.
In order to settle the nature of the transition between the perturbative and non perturbative regions we start by considering the radius of convergence of the Taylor expansion of around . It indeed measures the distance between the centre of the Taylor expansion and the nearest singularity in the complex plane. As the expansions of is an alternating series its possible singularity belongs to the real axis at the left side of the expansion centre. We can also argue on general grounds that the singularity should be at greater or equal to zero.
Indeed we expect at least a singularity of at , since measures the number of degrees of freedom of the subsystem . This number diverges exponentially with its size (Eq.(1.2), for example, presents an essential singularity at ).
For the sake of completeness one should remember that, when the system is in a product state (situation that does not arise close to a critical point) there is no singularity at , even in the thermodynamic limit. This happens, i.e., in the XY model described in (2.1) as one approaches the disorder line in the ordered phase . There the ground state becomes the superposition of two product states  and . Similarly, there is a class of quantum critical points in two spatial dimensions where a very simple form for has been proposed , namely , where and are the partition functions of the associated system with suitable boundary conditions which select the subsystem .
In the cases of highly entangled states we consider here, however, we do expect a singularity at in the thermodynamic limit and therefore the radius of convergence of the expansion (4.2) of centred around should not exceed 1. A precise value for it can be extracted using the ratio test. This states that if the limit exists, it coincides with the convergence radius . In a finite system, there is no true singularity and the radius of convergence of the series is infinite. However we expect we are still able to identify what would turn into a true singularity in the thermodynamic limit, by studying the sequence .
We anticipate that, notwithstanding there are clues to the presence of singularities in the thermodynamic limit, a careful analysis allows us to exclude it: the only possible singularities of lie on the imaginary axis of .
These ratio sequences show a sort of plateaux at a value in a certain range of before they eventually run away for . The point moves to the right within the size of the system, revealing longer and longer plateaux. This has been observed in all the models we have analysed. For instance, we may compute the required ratios from the Taylor coefficients of a block of a critical one dimensional Ising model, by using the coefficients plotted in Fig. 1 (a). These ratios are plotted in Fig. 3.
There we can see that the ratios oscillate for a while around a value smaller than 1 before running away. This suggests that a new singularity on the right side of the one we already expect at might exist in the thermodynamic limit.
This scenario is further supported by studying the Taylor expansions around other positive values of . If there were a true singularity on the real axis, the function
would converge to the coordinate location of the singularity in the thermodynamic limit and would not dependent on the value of . In a finite system we expect that should approach to a constant value independent of for a certain range of , and this range should increase with the size of the system. This behaviour can be observed, for instance, in Fig. 4 (a) where we plot for the block of the critical quantum Ising chain (the same already considered in Fig. 1 (a)), and in Fig.4 (b) where we plot the same quantity for the two-dimensional Ising model on a torus of side .
Still all the above signatures are not sufficient to conclude that there is a true singularity for real in the thermodynamic limit in .
In discrete systems of finite size partition functions are analytic with respect to their parameters because they are finite summations of positive terms. For the same reason they do not have zeros on the real axis. However there can exist zeros on the complex plane of certain control parameter, in our case, which are generally named Lee-Yang zeros . In general, the Lee-Yang zeros are apart from the real axis as long as the system size is finite. However, if a phase transition occurs at a critical value , the distribution of zeros becomes dense and accumulate on curves. The singularities of the free energy associated to this transition lie on the end points of these curves. As these end points move towards the real axis, eventually touching it at the location of the critical point .
In the case of besides the surviving zeros of the numerator there could also be poles due to zeros of the denominator. As an example, the exact expression of the reduced density matrix in Eq. (2.7) for the non-critical XY chain, presents a singularity at (the same as the singularity of (1.2)) but no zeros. It contains instead a huge number of double poles at for any pair of integers . This set is dense on the imaginary axis. The plot of is presented in Fig. 5 (a) and reflects this intricate structure. A similar landscape can be observed by considering the quantum two dimensional critical Ising model on a torus (see Fig. 5 (b)) where one can locate few zeros (see e. g. Figure 6) close to the region where the ratio test give hints of a singularity.
The important outcome of this study is that in all the models we have analysed both in and in , these Lee-Yang zeros, which are true singular points of the Rényi entropy (1.1) in the complex plane, do not become denser as the size of the system increases nor they approach the real axis (see Figure 7). In the thermodynamic limit then, the apparent singularity on the real axis identified by the ratio test of Taylor coefficients does not occur. The behaviour of flexes described by Eq.(5.1) is probably only a consequence of the presence of some structures in the complex plane (zeros and/or peaks) which do not evolve into a true singularity in the thermodynamic limit.
The lesson we learn from this analysis is twofold. On the one hand there is no true phase transition between the two regions where the Tsallis entropies show very different behaviour. The transition is only a smooth cross-over. On the other hand the behaviour of shows an intricate landscape in the region while it is very smooth in the complementary region . As a consequence a naive analytical continuation of the results obtained for the Tsallis entropies in the perturbative region is likely to fail to reproduce the correct entanglement entropy. We will provide an example of this fact in the following section.
As an aside, it is interesting to observe the striking resemblance between the location and shape of the zeros distribution of Fig. 7 and those found in certain disordered systems . In the latter case the replicas are coupled together no longer by the subsystem , but by the quenched average, which restores the translational invariance of the system. As a consequence the latter system develops in thermodynamic limit a true singularity associated to replica symmetry breaking.
6 The entanglement entropy mismatch
In order to test the consequences of the intricate structure of in the complex plane we have performed a numerical experiment. This consists in computing the entanglement entropy using the replica trick from the results obtained in a Monte Carlo simulation. For the 2D Ising model on a torus indeed there are still no analytical predictions on value of the parameter appearing in the proposed form of the scaling of the entanglement entropy. Some numerical results have been recently obtained in Ref. . Here we have computed for with the Monte Carlo technique described in Sect. 3 when is a half of the torus, the same geometry studied in Ref. . The results we obtain agree with the ones computed from the data of Ref.  provided we introduce a suitable conversion factor between the lattice spacings of these two different regularisation schemes. A good agreement is indeed obtained by multiplying the lengths measured in the Hamiltonian approach by a conversion factor (see Fig.8). This agreement is a good check of the validity of both methods. To our knowledge it is also the first direct comparison between two different non-perturbative techniques on entanglement entropy, one involving quantum mechanical tensor network simulations and the other one using statistical mechanics Monte Carlo integration. This is what we expect from universality, indeed, these two methods deal with very different systems (in one case the model studied is the 3D classical Ising model while in the other case it is the 2D Quantum Ising model). They show the same behaviour close to the critical point because they belong to the same universality class.
We can now try to extract a function from the data describing the behaviour of in the whole range of and considered . The function should depend on only few parameters and is chosen in such a way to agree, in the limit, with the expected general form of the critical entanglement entropy in 2+1 dimensions [23, 10]. A sufficiently good fit (a reduced for 17 degrees of freedom) is given by the three-parameter function
The ‘obvious’ extension to real , is performed by assuming that this functional form should also be valid for real . The entanglement entropy then is given by
Due to our choice of the parametrisation (6.1) the entanglement entropy has the expected functional form The numerical coefficients we extract however, apart from the ’area term’ , do not coincide with the ones extracted from the direct study of the entanglement entropy computed in Ref. . If we call the parameter extracted in Ref.  , taking into account the conversion factor , we would expect from the values reported in (6.2) that
On the other hand their direct numerical determination  provides different results
We attribute the mismatch of the entanglement entropy to the the different behaviour of the the regions and .
A further confirmation of this can be obtained by computing, in the Hamiltonian approach, for some real (or complex) in the non-perturbative region. If we add to the fitted data some of these points the quality of the fit to (6.1) rapidly deteriorates.
One can check that a similar phenomenon holds also for models. Indeed, considering again the critical XY system in which the subsystem is composed by a block of adjacent spins, can be evaluated with the method of  and the data with can accurately fitted with Eq.(1.2) ignoring the dependence of the non-universal coefficient which is known to be neglegible in such a range . If now we assume that this functional form with is also valid for real and try to evaluate this quantity near we obtain incorrect results for the entanglement entropy. As an example the Taylor coefficients of this fitting function are plotted in Fig. 9 and should be compared with the same coefficients directly evaluated on (see Fig.1 (a) ). They show a completely different shape. Again this mismatch has to be attributed to having applied to the fitting function (1.2) with , which accurately describes the behaviour of in the perturbative region (i.e. ), the ‘obvious’ extension to real close to which lies in the non-perturbative region.
In this work we analysed the trace of the th power of the reduced density matrix associated to a subsystem in some Ising-like quantum models in both one and two spatial dimensions near and at a quantum critical point at zero temperature. We unveiled that, depending on the value of , shows two very different behaviours. If is larger than a threshold value -located in the interval - its behaviour is very smooth. In particular the alternating coefficients of its Taylor expansions, are monotone decreasing. In this regime the error involved in a truncation of the series does not exceed the first excluded term and can be always chosen to be very small. We label this region ‘perturbative’ region.
In the range , on the other hand, the ’s are no longer monotone decreasing function of but they present a characteristic peak (see Fig.1 ) whose height increases with the size of the block considered. Any truncated Taylor expansion, in this region, provides an increasingly poor approximation of as the size of increases. This is the region that we label ‘non-perturbative’.
The transition between these two regions turns out to be a smooth cross-over. The dependence of the threshold on the size of the subsystem obeys a power law described in Eq. (4.6). We determined the critical exponents ruling this power law for the Ising model in both one and two dimensions. It would be interesting to investigate whether this cross-over is somehow related to the phase transition close to n=2 recently observed in  near .
The above scenario has important consequences when considering the analytical continuation involved in the computation of the entanglement entropy with the replica trick. This continuation is possible in one and two dimensions, as there is no phase transition between the perturbative region (where both analytical results from QFT and numerical results from Monte Carlo algorithms are available) and the non-perturbative region (where the entanglement entropy is located). However its form is not the naive expression obtained by promoting the simple expression obtained in the perturbative region for integer to general , as we have explicitly shown in Section 6 .
From our analysis, the Rényi and Tsallis entropies seem to emerge as a privileged measures of entanglement in QFT. Although they lack the nice theoretical information properties of the entanglement entropy , their analytic properties as functions of the parameter seem to contain some new information on the entanglement. This has been considered recently in the context of quantum criticality [29, 30, 10] and topological order .
LT acknowledges the financial support from the QuSim group of the University of Queensland, the stimulating discussions with P. Corboz and the enlightening comments on entanglement measures of Prof. G. Vidal.
-  L.Amico, R.Fazio, A. Osterloh and V. Vedral, Rev. Mod. Phys. 80, 517 (2008), arXiv quant-ph/0703044; J. Cardy, Eur. Phys. J. B 64 , 321 (2008), arXiv 0708.2978; J. Eisert, M. Cramer and M.B. Plenio, arXiv 0808.3773.
-  A. Rényi, “ Probability Teory”, North-Holland, Amsterdam,1970.
-  C. Tsallis, J. Stat. Phys.52 (1988) 479.
-  G. Evenbly, G. Vidal Phys. Rev. Lett. 102, 180406 (2009) arXiv:0811.0879.
-  F. Verstraete, J.I. Cirac, V. Murg Adv. Phys. 57,143 (2008) arXiv:0907.2796.
-  C. Holzhey, F. Larsen, F. Wilczek, Nucl.Phys. B424 (1994) 443 arXiv:hep-th/9403108
-  P. Calabrese, J. Cardy J.Stat.Mech. 0406 (2004) P002 arXiv:hep-th/0405152. P. Calabrese, J. Cardy arXiv:0905.4013
-  J. L. Cardy, O. A. Castro-Alvaredo, and B. Doyon J. Stat. Phys. 130 1 (2008) arXiv:0706.3384.
-  JJ.L. van Hemmen and R.G. Palmer, J. Phys. A: Math. Gen. 12, 563 (1979).
-  M.A. Metlitski, C.A. Fuertes and S. Sachdev, arXiv:0904.4477.
-  G. Vidal, J. I. Latorre, E. Rico, A. Kitaev, Phys. Rev. Lett. 90, (2003) 227902 arXiv:quant-ph/0211074.
-  A. R. Its, V. E. Korepin arXiv:0906.4511.
-  L. Tagliacozzo, G. Evenbly and G. Vidal, arXiv:0903.5017.
-  M. Caraglio and F. Gliozzi, JHEP 11,076 (2008) arXiv:0808.4094.
-  H. Rieger, N. Kawashima, Europ. Phys. J. B 9, 233 (1999) arXiv:cond-mat/9802104 . H.W.J. Blöte and Y. Deng, Phys. Rev. E 66, 066110 (2002).
-  B.-Q.Jin and V.E. Korepin, J. Stat. Phys. 116,79 (2004) arXiv:quant-ph/0304108.
-  A.R. Its, B.-Q.Jin and V.E. Korepin, Journal Phys. A 38, 2975 (2005) arXiv:quant-ph/0409027.
-  V. Popkov and M. Salerno, Phys. Rev. A 71, 012301 (2005) arXiv:quant-ph/0404026.
-  T. Nishino, J. Physs. Soc. Jpn. 74, 3598 (1995) arXiv:cond-mat/9508111.
-  I. Peschel, J. Stat. Mech. P12005 (2004) ArXiv:cond-mat/041041.
-  F. Franchini, A.R. Its and V.E. Korepin, J. Phys. A 41 025302 (2008) arXiv:quant-ph/0707.2534.
-  P. Calabrese and A. Lefevre, Phys. Rev. A 78, 032329 (2008) arXiv:0806.3059.
-  H. Casini and M. Huerta, Nucl.Phys.B 764 183 (2007) arXiv:hep-th/0606256.
-  I. Peschel and V.J.Emery, Z. Physik B 43 241 (1981).
-  E. Fradkin and J.E. Moore, Phys. Rev. Lett.97,050404 (2006), arXiv:cond-mat/0605683.
-  C.N. Yang and T.D. Lee, Phys. Rev. 87, 404 (1952); T.T. Lee and C.N. Yang, Phys. Rev. 87, 410 (1952).
-  K. Ogure and K. Kabashima, J. Stat. Mech. P05011 (2009), arXiv:0812.4655.
-  C. H. Bennett, H. J. Bernstein, S. Popescu, and B. Schumacher Phys. Rev. A 53, 2046 - 2052 (1996) arXiv:quant-ph/9511030.
-  J.-M. Stéphan, S. Furukawa, G. Misguich, V. Pasquier arXiv:0906.1153.
-  E. Fradkin arXiv:0906.1569.
-  S. T. Flammia, A. Hamma, T. L. Hughes, X.-G. Wen arXiv:0909.3305.