Entanglement at the quantum phase transition in a harmonic lattice

The entanglement properties of the phase transition in a two dimensional harmonic lattice, similar to the one observed in recent ion trap experiments, are discussed both, for finite number of particles and thermodynamical limit. We show that for the ground state at the critical value of the trapping potential two entanglement measures, the negativity between two neighbouring sites and the block entropy for blocks of size 1, 2 and 3, change abruptly. Entanglement thus indicates quantum phase transitions in general; not only in the finite dimensional case considered in [Phys. Rev. Lett. {\bf 93}, 250404 (2004)]. Finally, we consider the thermal state and compare its exact entanglement with a temperature entanglement witness introduced in [Phys. Rev. A {\bf 77} 062102 (2008)].


Introduction
Coupled harmonic chains with short-and long-range (LR) interactions are ubiquitous in science and engineering. Their application to calculate the phononic heat capacity by Einstein [3] marks the birth of solid state physics. Beyond physics, harmonic chains feature in chemistry and biology, where they are used to model the behaviour of macro-molecules, such as DNA [4] and cell membranes [5]. In the last decade, harmonic systems have been revised using techniques developed in quantum information science to study correlation properties in the quantum regime and particularly at small temperatures [6]- [9]. Thermodynamics has been very successful in characterizing 'standard' phase transitions that occur at finite temperature when a macroscopic parameter, such as pressure, is changed [10]. Quantum phase transitions (QPTs) appear at zero temperature [11] and are due to the change of an external parameter, such as the trapping potential of an ion trap. These transitions are driven by quantum fluctuations and have been linked to entanglement for the case of finite-dimensional systems [1,12].
In this paper we study a QPT in a continuous variable system: a system of trapped ions that we model as a harmonic lattice [13,14]. The ions interact via an LR Coulomb repulsion and are trapped by two external potentials; see figure 1. They align in a linear configuration for sufficiently large transversal trapping potential ν t ; however, when ν t is decreased the system undergoes a phase transition and the new equilibrium state forms a zig-zag configuration. This model is motivated by ion trap experiments [15]- [17], where such a QPT occurs [18,19]. Recent analytical studies of the transition using Landau theory [20] allowed one to determine the system's classical behaviour at the transition point. Moreover, the numerical treatment of the quantized system of a few ions promised the possibility of simulating linear and nonlinear Klein-Gordon fields on a lattice [21]. However, a comprehensive analytical study of the quantized system has so far been lacking due to the complexity of the system.
Here we model the ion trap scenario as a lattice of harmonically coupled oscillators and present a first quantitative characterization of the entanglement inherent in both ground state configurations of the ions. For finite-dimensional systems, QPTs of first (second) order are characterized by a discontinuity (a discontinuity in or divergence of the first derivative) of the negativity [1]. We also show that in the continuous variable system considered here the structure Sketch of the harmonic lattice under consideration: each site is trapped by two external potentials, ν in the x-direction and ν t in the y-direction. For clarity only NN coupling is indicated; however, all ions interact via an LR potential, which we approximate harmonically. This could be, for instance, a Coulomb potential such as is common in ion experiments. The equilibrium distances between sites are assumed to be equidistant, with lattice constant a in the x-direction and b along the y-direction. In I, the transverse trapping potential ν t is larger than the critical value and the ions arrange linearly. As displayed in II, decreasing the transverse trapping potential below the critical value ν t,crit leads to a QPT causing the ions to move outwards and form a two-dimensional zig-zag structure. of entanglement, measured by the negativity and the von Neumann entropy, changes abruptly at the critical point and indicates the occurrence of a QPT, in a manner similar to the way classical correlations indicate standard phase transitions. The first derivative of the negativity between two neighbouring ions has a finite discontinuity and the von Neumann entropy of contiguous blocks of a single ion, two and three ions all show a divergence in the first derivative. The LR nature of the Coulomb interaction leads to an increase of the block entropy with increasing block size. This is in contrast to models with only nearest neighbour (NN) interaction where 'area laws' apply [22] and entanglement does not increase with block size, i.e. volume, as long as the surface of the block is constant. However, our results show that the increase in block entropy with block size is quite small due to the fast decline of the Coulomb potential.
The quantum fluctuations that cause the QPT are most dominant at zero temperature. However, as experiments are performed at small but finite temperatures, it is important to know how temperature affects these fluctuations [23]. In the final part of our paper, we discuss thermal states and find that the sharpness of the QPT, indicated by the entanglement, fades out with increasing temperature. Another macroscopic consequence of quantum fluctuations is the lowering of the energy of the system [2]. We compute up to which temperature the thermal state 4 has a lower energy than any separable state. This temperature witness could be implemented in an ion trap experiment by measuring the average energy (i.e. mean excitation). Although the model is motivated by ion traps, realizations include the vibrational motion of molecules [4] (nuclei in the electronic potential) and optical lattices [24].

The model
We consider the Hamiltonian (see figure 1) where are the Coulomb potentials of the sites with r j = (x j ,ỹ j ) absolute coordinates of the sites andx j,0 are the equilibrium positions in the x-direction. N is the number of particles, Q is the charge of the ions and ν (ν t ) are the trapping potentials in the x (y)-direction. We assume periodic boundary conditions, r j = r j+N .
To calculate the entanglement measures, we approximate the Coulomb potential to second order and expand about the equilibrium positions. The key step is then to diagonalize the Hamiltonian into a set of uncoupled modes, the lattice vibrations, with which analytic expressions for the measures can be obtained. Similar to the classical calculation [20], we use a simplified model with equidistant equilibrium position in the x-direction, spaced by the lattice constant a. Such a condition can be realized for the central ions of a long ion chain inside a linear Paul trap [25] or for ions confined in a ring of large radius [18,26].
For large trapping potential the sites are arranged on a single line, i.e.x j,0 = a j andỹ j,0 = 0, while for small enough ν t the equilibrium positions becomex j,0 = a j andỹ j,0 = (−1) j b/2 and a two-dimensional zig-zag configuration is formed. The equation determining deviation b in the y-direction is obtained by summing the linear terms over all sites and requiring it to vanish, 1 2 where τ = k − j numbers the neighbours of each site. The harmonically approximated Hamiltonian becomes with x j =x j −x j,0 and y j =ỹ j −ỹ j,0 the deviations from equilibrium. Furthermore, the d x,y,x y denote second order Taylor coefficients of the Coulomb potential which are, for the linear and zig-zag configuration,

Calculation of entanglement measures
We are interested in the behaviour of the entanglement between sites in the chain for varying transverse trapping potential, ν t , and particularly at the point of criticality, ν t,crit . We calculate the negativity, E N , for two modes regardless of all others, for instance, the x degrees of freedom of two (neighbouring) sites. To measure the correlation of one mode with all other modes, for instance the entanglement between the y degrees of freedom of a single site or block of sites with all the other degrees of freedom in the chain, the von Neumann entropy, S V , is used. Other entanglement measures are available, such as the entanglement of formation [27]. Yet they are very hard to calculate in this continuous variable scenario and we will be content with the two measures as stated. To see the effect of the LR interaction, we compare the full LR Coulomb Hamiltonian with a cut-off version in which only NNs interact and the interaction with more distant neighbours is set to zero.

Diagonalization of the Hamiltonian
For the linear configuration, the Hamiltonian decouples into x and y parts. A discrete Fourier transformation for x (similar in the y-direction) of the form maps the space coordinates of the sites into diagonal modes, lattice vibrations or phonons. The diagonal Hamiltonian is then H = N l=1h ω xl n xl + 1 2 + N l=1h ω yl n yl + 1 2 , wheren xl = (P † xl P xl /2h mω xl ) + (mω xl X † l X l /2h) − 1 2 and similarlyn yl are the number operators in the xand y-direction for mode l. The frequencies are where C = 4Q 2 /ma 3 . The asymmetry of the system is reflected in the dispersion relations. ω xl is always real, whereas ω yl becomes complex for small values of transverse trapping potential ν t . This is where the quantum phase transition occurs. The critical value ν t,crit for LR interaction is ν t,crit ≈ √ 0.6C whereas for NN interaction it is ν t,crit = √ 0.5C. To diagonalize the Hamiltonian in the emerging zig-zag configuration, where previously independent phonons now couple in the x-y-direction, see equation (3), we need to amend the transformation in the y-direction (5) and transform with 6 The additional factor e iπ j = (−1) j compensates the alternating sign of the x-y coupling in the Hamiltonian. Expressed with a coupling matrix the Hamiltonian can now be written as with coefficientsω u,l = ν 2 For each l we need the symplectic transformation [28], denoted by S l , that diagonalizes the matrix M l , where x yl (14) are the symplectic eigenvalues of M l . These transformations are given by where and using the number and similarly for w, we find the fully diagonalized Hamiltonian for the zig-zag configuration H = l h ω v,l n v,l + 1 2 +hω w,l n w,l + 1 2 .

Two-site negativity and von Neumann entropy
The calculation of the two-site negativity requires the evaluation of the covariance matrix of the partially transposed state of the two sites [13,14]. In two dimensions this is an 8 × 8 matrix of which the symplectic eigenvalues have to be found. However, in the linear configuration of the ions, x and y degrees of freedom completely decouple and it is sufficient to consider two 4 × 4 covariance matrices independently. In contrast, the zig-zag configuration contains x y-coupling terms and it is a priori necessary to consider the full 8 × 8 matrix. As a result, no entanglement occurs between the xand y-direction in the linear configuration while the zig-zag configuration can sustain x y-entanglement. However, we found that all expectation values coupling the xand y-direction, i.e. x i y j , p x i p y j etc, vanish also in the zig-zag configuration.
To characterize the entanglement between two sites, two sets of each two conditions, as used in [2], emerge. For each site j these separability conditions are and similarly for the y-direction. The expectation values needed here can be calculated using the transformation rules into the diagonal modes, equation (7) in the linear chain and equations (9) and (16) in the zig-zag configuration. If one of the inequalities is violated, then entanglement exists between the jth site and its τ th neighbour and the negativity measures their degree of entanglement. The two criteria S 1 and S 2 witness two types of entanglement. For example, the Einstein-Podolsky-Rosen pair originally considered in [30] shows violation for S 2 but not for S 1 .
The von Neumann entropy of a single site j in either x or y dimension is obtained following [31] with the formula where r j = x 2 j p 2 x j , and similarly for the y-direction, is the symplectic eigenvalue of its reduced state. To evaluate the entropy of a block of n neighbouring sites, the block entropy is then simply where the sum is taken over all n symplectic eigenvalues r j in either x or y dimension within the block.

Thermodynamical limit (N → ∞)
To obtain the negativity in the x-direction at zero temperature in the thermodynamical limit in the linear configuration, we evaluate equation (17) in the ground state leading to and similarly for the entanglement in the y-direction. For the zig-zag configuration transformation (16) gives a more complicated expression for the S 1,2 criteria. The formula is long and not enlightening; therefore we omit it. In the thermodynamical limit we substitute πk/N = α, πl/N = β, and replace the sums with integrals, These integrals can be evaluated numerically and the plots are shown in figure 2 (red-solid line). For one special set of parameters the above expression can be easily calculated analytically, namely for the y-direction in the linear configuration at the critical point ν 2 t,crit = C/2 for NN interaction. The S 1 criterion between NNs becomes which leads to a value of the negativity of E N ≈ 0.308. In a similar fashion, the single-site von Neumann entropy for both xand y-directions can be evaluated in the thermodynamical limit. Here we show the calculation for the ydirection. The symplectic eigenvalue for a single site in the y-direction is r j = y 2 j p 2 y j = N k,l=1 (1/N 2 ) Y 2 l P 2 y l . Again we substitute sums with integrals and πk/N = α, which gives an integral expression for the symplectic eigenvalue At the critical point ν 2 t,crit = C/2 and for NN interaction, this can be simplified to The first integral diverges and hence the symplectic eigenvalue and the von Neumann entropy diverge at the QPT. This is because phonons in the x-direction are independent of trapping in the y-direction. At the critical point both negativity and entropy are not differentiable. Decreasing the trapping potential beyond the critical value, where the zig-zag configuration is formed, the x-entanglement is reduced due to the emerging x-y coupling. Because each site in (b) couples to several different neighbours, the entanglement between two distinct sites disappears.

Behaviour of entanglement at zero temperature
In the lower plots (c and d) both entanglement measures for the y-entanglement increase with decreasing trapping potential ν t in the linear configuration. The even-numbered ions oscillate exactly out of phase with the odd numbered ions due to the repulsive Coulomb potential. The smaller the ν t , the larger these quantum fluctuations around the equilibrium position. At the critical point the fluctuations become strong enough for causing the ions to move outward. The entropy diverges and the negativity reaches its maximal value of E N ≈ 0.308 where it is not differentiable.
For NN coupling (a and c), the negativity and entropy qualitatively show the same behaviour. This can be understood easily as any entanglement of a site with the rest of the chain, measured by the von Neumann entropy, is created by coupling with only the NNs. For LR interaction (b and d), where a single site couples to all other sites, there are significant differences between the two entanglement measures. While the negativity of two NN sites vanishes after a threshold value of ν t , each single site remains entangled with the rest of the chain, as seen by the positive value of the von Neumann entropy.
Both entanglement measures are functions of the eigenvalues of the Hamiltonian. Therefore abrupt changes in entanglement can signal the non-analyticity of ground state energy, which is associated with QPTs. In [1] it was shown for finite-dimensional systems that (under certain conditions) a discontinuity in or divergence of the first derivative of the negativity is both necessary and sufficient to signal a QPT. It seems intuitive that a similar characterization also holds for continuous variable systems. Here the divergence of entropy and finite discontinuity of the first derivative of negativity shows a QPT of second order. After the phase transition, the violated negativity criterion switches from S 1 to S 2 . Lowering the trapping potential further leads to another 'critical' point c y at which the negativity becomes zero while the entropy reaches its minimal (nonzero) value.
The additional critical point c y (and c x in (b)) only appears when the interaction is harmonic as is the case in our second order approximation of the Coulomb model. The point is due to a sign change of the second order coefficient d y τ (d x τ ) at c y (c x ). As a consequence, the interaction switches from repulsive to attractive. When anharmonic terms are taken into account, as in the numerical treatment in [21], these points vanish. An intuitive way of understanding the switching between S 1 and S 2 is as follows, see figure 3. Decreasing the trapping potential changes the relative strengths of the inner and outer potential for motion in the y-direction. This leads to the change in phase in the relative momenta of neighbouring sites. In one configuration the relative potential favours momenta in the opposite directions, whereas the other configuration favours motion in the same direction. This is reflected in the change from S 1 to S 2 . Although the negativity vanishes at these points, a single site is still entangled with the rest of the chain.
Block entropy measures how much entanglement exists between a block of sites of the lattice and the rest. For NN interaction models, there exist scaling laws showing that the amount of entanglement scales with the boundary area and not the volume of the reduced state [22]. For a translational invariant chain with NN interaction, there exists even a computable analytical result for the negativity of a bisected harmonic chain [6]. Here we investigate block entropy for   blocks up to three sites in the LR Coulomb lattice. Due to the complexity of LR interactions, few results are known so far. Inset (c) illustrates the increase of correlation across the block boundary with increasing block size. Correlations stretch to NNs and NNNs, as shown in figure 4, while third neighbour entanglement is negligible. Increasing the block size from one to two, there are twice as many NNN connections across the boundary and hence the block entropy for two sites is expected to increase. This intuition is confirmed in (a) and (b), showing block entropy in the xand y-direction, respectively. However, it can be seen that while the entropy increases slightly with number of sites, no qualitative difference can be observed. This is due to the fact that the additional, LR entanglement is much weaker because the Coulomb potential falls off quickly with 1/τ 3 , where τ is the distance of sites. In y dimension (b), already the single site entropy is a good approximation for larger blocks. This is because transversal NNN entanglement turns out to be very small. In the x-direction (a), there are significant differences for growing block sites.
Here NNN entanglement cannot be neglected. Although the number of connections between two and three sites is the same, block three entropy is still higher. This might be due to multipartite entanglement. Larger block sites are difficult to evaluate, as the symplectic eigenvalues of the covariance matrix become very complicated.

Witnessing entanglement at finite temperature
One consequence of the entanglement of sites is a lowering in the energy of the system [2]. This can be seen by assuming that the thermal state of the system is separable, i.e. decoupled, between modes. Then an effective, single site Hamiltonian can be obtained by removing all second order couplings between the different sites, i.e.
The total Hamiltonian becomes the sum of the single site Hamiltonian where Then the thermal state takes the form with β = 1/k B T the inverse temperature. Using the transformation, equations (7), (9) and (16), H eff can be fully diagonalized. The internal energy U = H eff for any separable state is bounded from below by zero point fluctuations, i.e.

H eff sep
Nh and any state having a smaller energy must be entangled between individual sites. We now want to see how the ground state situation is modified at nonzero temperature. As the energy of the thermal state, i.e. the mean excitation of phonons, increases with temperature, there exists a critical temperature, T c , for each value of the trapping potential at which the thermal state matches the energy bound. The negativity between the y degrees of freedom of two nearest neighbours for the NN Hamiltonian is evaluated numerically and plotted in figure 5, where the critical temperature for full separability is also indicated as a red line. When the trapping potential is lowered, the negativity increases until its maximal value at the critical ν t,crit . Further lowering leads to a decrease of the negativity until it vanishes at point c y . However, when further decreasing ν t the negativity increases again, yet the two entanglement criteria S 1,2 are now switched. As expected, increasing the temperature leads to smaller values of negativity and smoothens the entanglement measures to make it differentiable at the critical point. For large ν t the negativity is small, but remains finite until relatively high temperature. The sharp peak at the QPT remains almost constant for finite T and decreases rapidly. Thermal states within the red outlined area have a smaller energy than any separable state and their entanglement is therefore detected by the energy witness. Remarkably, states at c y are entangled for temperatures up to T c [1/ν t = c y ] = 0.12[ν t ]h/2k B , even though there is no NN entanglement in the y-direction. This is because there is still NN entanglement in the x-direction and possibly also multipartite entanglement in the chain. . The red line indicates critical temperature, obtained with the energy witness argument, below which entanglement of some form must be present in the chain.

Conclusions
In this paper, we revised the classical phase transition in an LR harmonic chain [20] using a fully quantized model. Two measures of entanglement display critical behaviour: the von Neumann entropy of a single site and blocks of two and three sites diverge at the critical point while the negativity is not differentiable. Thus also in this continuous variable system, entanglement indicates a QPT, as previously shown for discrete systems [1]. The negativity depends only on single site and NN correlations; the single site von Neumann entropy only depends on single site measurements. Our calculation shows that even this local, single site function is able to detect a global change in configuration. This implies that instead of examining two point correlation functions, one can alternatively consider the entanglement measures stated. This will be advantageous in experimental situations when the number of different measurement procedures is best kept as low as possible. We are aware that for the moment ion traps cannot yet perform the required measurements of e.g. single site variance of the space and momentum operator, but this is a technical, not a fundamental, problem. Furthermore, our results confirm that this phase transition is of second order, as indicated in [20]. At finite temperature, the negativity still displays critical behaviour, as seen in figure 5; however, the non-differentiable cusp fades out quickly with increasing thermal noise. Tuning across a QPT provides a means of changing the amount and structure of the continuous variable. Experiments with ion traps are ideally suited to study QPTs with great precision.