Entanglement in Lifshitz-type quantum field theories

We study different aspects of quantum entanglement and its measures, including entanglement entropy in the vacuum state of a certain Lifshitz free scalar theory. We present simple intuitive arguments based on “non-local” effects of this theory that the scaling of entanglement entropy depends on the dynamical exponent as a characteristic parameter of the theory. The scaling is such that in the massless theory for small entangling regions it leads to area law in the Lorentzian limit and volume law in the z → ∞ limit. We present strong numerical evidences in (1+1) and (2+1)-dimensions in support of this behavior. In (2 + 1)-dimensions we also study some shape dependent aspects of entanglement. We argue that in the massless limit corner contributions are no more additive for large enough dynamical exponent due to non-local effects of Lifshitz theories. We also comment on possible holographic duals of such theories based on the sign of tripartite information.


Introduction
Studying the physics of non local correlations due to the quantum entanglement in quantum many-body systems, quantum field theories and especially holographic field theories has gained increasing attention during the last decade [1][2][3][4]. Furthermore, in order to quantify entanglement and have a deeper understanding of it, different measures has been studied so far such as entanglement entropy (EE), mutual information and logarithmic negativity [5][6][7][8]. In particular EE is a good measure of entanglement for pure stets, since it behaves as an entanglement monotone decreasing under LOCC [9].
The recipe for computing EE is simple and straightforward. Considering the simplest setup where the physical system consists of two subsystems such that H tot. = H A ⊗ HĀ. Integrating out the degrees of freedom which live in the second subsystem, i.e.,Ā, one can find a reduced density matrix for A which we denote by ρ A . EE is given by the corresponding von Nuemann entropy for ρ A . In contrast to the thermal entropy, the entanglement entropy is nonvanishing at zero temperature and thus it is a good probe for studying properties of quantum phase transitions. Also it captures some characteristic properties of the underlying theory which the best known examples are even dimensional conformal field theories (see [10,11] for reviews).
At a quantum critical point, a physical system typically exhibits a Lifshitz scaling symmetry as follows [12][13][14] t → λ z t, x → λ x, (1.1) where z denotes the dynamical critical exponent. The above inhomogeneous scaling in general breaks down the beloved Lorentz invariance. Merely in the very special case of z = 1,

JHEP07(2017)120
the relativistic dynamics is recovered. A huge amount of efforts has been made to understand the entanglement structure in relativistic QFTs, which has resulted in many interesting features and of course many open questions. This constructs the very most literature of quantum entanglement in the context of QFTs. Although there have been some attempts so far to analyse the entanglement structure in the non-relativistic cases with z = 1 (see [15][16][17] on the QFT side and [18][19][20][21][22][23] in the context of gauge/gravity duality 1 ), still there are several open questions about the role of z in the entanglement structure of such theories.
Here one of our main motivations is to construct a more stringent analysis of the nature of quantum correlations at Lifshitz critical points. An interesting tool to do so is studying quantum entanglement in theories with Lifshitz symmetry. We will consider the simplest example of such QFTs, which is a free scalar theory with generic dynamical exponent z. This family of QFTs represents several interesting features which is not of our interest here (see e.g. [28] for a review). A similar analysis has been previously done for the same theory but for 0 < z < 1 in [29][30][31][32]. We show that the behavior of z > 1 is totally different from what was found for z < 1. 2 The first step to study the entanglement structure of these theories is to focus on wellknown measures such as entanglement entropy. This is often a formidable task to carry out analytically in the context of QFTs. In principle there are different strategies including, replica trick [2], numerical approaches or even 'generalized' holographic prescriptions [33] which one can employ them to overcome this problem and gain some information about the entanglement structure of theories with Lifshitz symmetry.
In this paper we will tackle the problem numerically using the so called correlator method [34,35]. This method is based on a semi-analytic approach leading to the reduced density matrix (as well as its n-th power) corresponding to a given subregion. The correlator method which is applicable to quadratic field theories is based on reducing all correlation functions to two-point functions restricted to a given subregion. This method which we will shortly review in the next section makes it possible to calculate entanglement and Renyi entropies in free field theories. We will employ this method to study entanglement measures in the ground state of a free scalar theory with Lifshitz symmetry.
We will focus on free massless and massive scalar theory with dynamical exponent z ≥ 1. In order to perform numerical calculations we put the theory on a square lattice. The main difference between theories with Lifshitz symmetry and their Loretzian counterpart, i.e. the same theory with the dynamical exponent z = 1, is that the number of lattice sites which are 'correlated' together due to the spatial part of the kinetic term is z-dependent.
To be more precise, for a positive integer z, the kinetic term has a spatial derivative of order 2z, thus it correlates (2z + 1) lattice points together. This fact revises the well-known intuitive understanding of the area law of entanglement in local field theories which is based on the short range correlation due to the kinetic term. In figure 1 we have illustrated the intuitive counterpart in this case that which the area law transits to a volume law while the dynamical exponent increases. 1 QFTs with Lifshitz symmetry has been widely studied in the context of gauge/gravity duality, see e.g. [24][25][26][27]. 2 We would like to thank Tadashi Takayanagi for bringing our attention to these references during the final steps of this work.

JHEP07(2017)120
. . Figure 1. The nearest points inside a disk entangling region to its boundary are shown in green. Those points which are correlated to the green ones due to the discrete kinetic term of this theory are shown with a red shadow. The left panel belongs to the Lorentzian case z = 1. This shows how we intuitively understand the area law of entanglement for regions which their characteristic length, say , is much bigger than the lattice spacing (the inverse of the UV cut-off). Moving from the left panel to the right, the dynamical exponent is increasing. One can see that for a disk with a radius ∼ 3.5 in units of lattice spacing, for z > 6, since all points inside the entangling region are correlated with the green ones, we expect entanglement measures to scale with the volume (instead of the area) of the entangling region. We will show numerically that this is the correct scaling of entanglement entropy for large enough z in the following of this paper. Although this figure belongs to (2 + 1)-dimensions, its horizontal (vertical) slices describe what happens in (1 + 1)-dimensions and it can be generalized to higher dimensions straightforwardly.
The rest of this paper is organized as follows: in section 2 we briefly introduce our QFT setup and correlator method which we use in our analysis. We report our results for different entanglement measures of our Lifshitz theory of interest in the massless regime in section 3 and show that EE obeys a volume law for large z. We also study the shape dependence of EE in 2+1 dimensions and report our results for different entangling regions such as a disk and square and present strong evidences for emergent volume law behavior of EE at large z. We also show that corner contributions to EE get away of being local effects in such theories for z > 1. In section 3.2 we continue our study by considering a massive scalar theory in different dimensions. The last section is devoted to the conclusions and several possible interesting directions to investigate in future works.
Note added. Reference [36] which appeared after our paper has some overlap with our study in (1+1)-dimensions. The results are compatible where overlap exists.

QFTs with Lifshitz symmetry
The starting point of our analysis is the following nonrelativistic action for a free massive scalar field in (d + 1)-dimensions [28]

1)
3 Note that the most general theory with such a symmetry may also include terms with lower order spatial derivatives, say φ (∂ · ∂) ν φ with ν being an integer and ν < z. These terms have generic coefficients in the action which (for simplicity) we have set all of them to vanish in our analysis (see [28] for detailed analysis of such theories).

JHEP07(2017)120
which in the massless case is invariant under the scaling introduced in eq. (1.1). According to this action the mass dimensions are given as follows where z should be an integer parameter. In the following discussion we always consider z > 1 case. In this theory in the massless case one can easily find the scaling behavior of the ground state two point correlator to be [32] This behavior shows that the power law growth of correlator as the dynamical exponent increases. This is in agreement with what we have explained intuitively in figure 1.
Using the momentum density conjugate to φ, i.e., π = ∂L ∂φ , the corresponding expression for the Hamiltonian density can be written In the following we concentrate on case d = 1 having in mind that the generalization to arbitrary d is straightforward. In order to push our calculations we introduce a UV cut-off and thus put the theory on a lattice. The Hamiltonian density on a square lattice is given by where we consider a lattice with N sites and without loss of generality we set the lattice spacing equal to unity. Note that φ n and π n satisfy standard commutation relations, i.e., [φ n , φ m ] = [π n , π m ] = 0 and [φ n , π m ] = iδ nm . In the following sections we will consider periodic and Dirichlet boundary conditions which correspond to φ N = φ 0 (π N = π 0 ) and φ N = φ 0 = 0(π N = π 0 = 0) respectively. In the case of periodic boundary conditions, using the translation invariance one can diagonalize the Hamiltonian given in eq. (2.5) with the following Fourier transformations Expanding these new canonical variables in terms of creation and annihilation operators, i.e.,φ the Hamiltonian density becomes

JHEP07(2017)120
where (2.9) We have used the following sum rule (2.10) to derive the dispersion relation given in eq. (2.9). Notice that although in eq. (2.1) we consider z as an integer parameter, actually eq. (2.9) shows the exact analytic continuation to non integer dynamic critical exponents. Employing the above definitions of dynamical fields in terms of creation and annihilation operators the vacuum correlators are given as follows where for the infinite lattice limit, i.e., N → ∞, becomes where Ω(x) 2 = m 2z + (2 sin πx) 2z . Of course the correlators only depend on (m − n) due to the translation invariance of the system. 4 Note that according to the above expressions, in the case of periodic boundary condition the massless limit, i.e., m → 0, is not well-defined due to the existence of a zero mode. This makes us consider a non-zero mass in order to regularize the corresponding divergences. One way to get rid of the zero mode is to break the translational symmetry of the system. As an example one can consider Dirichlet boundary condition on the borders of the total system. Since there is no translation invariance in this system, one should consider a Fourier sine transformation for the fields which gives [38], The analysis then proceeds in analogy with the case of periodic boundary condition. Evaluating the corresponding two point functions, we find that In this case the zero mode is no more excited and there is no divergence in eq. (2.14) while m → 0. In this case we will consider the massless limit. 5 4 See [37] for similar analysis of Lorentzian free scalar theory. 5 An explicit analysis of zero modes in EE in relativistic free scalar theory has been recently carried out in [39].

Review of correlator method
As we mentioned in the introduction section, in order to find the entanglement entropy, one should construct the reduced density matrix corresponding to the subregion of interest. Since we are interested in a free (quadratic) theory, we use the correlator method first introduced in [34]. 6 This method is based on a very simple idea: if we are with the expectation value of any operator restricted to region A, which we denote by O A , in principle this would fix the reduced density matrix corresponding to this subregion uniquely through In the context of quadratic theories, thanks to Wick's theorem one can reduce the right hand side of eq. (2.15) to all possible 2-point functions of the theory. In this way it is not hard to show that the eigenvalues of the reduced density matrix which has the following form where ν k 's are the eigenvalues of C = √ X.P which itself is given by the following 2-point functions Using the definition of entanglement and Renyi entropies it is a very simple task to show that where n A is the number of lattice points enclosed by the entangling region A. The method just outlined is very general and can be employed to any quadratic QFT with different matter content such as fermions in any space-time dimensions. Also note that using this method one can find other entanglement measures, e.g., mutual and tripartite information as we will do in the subsequent sections. restricts us to impose Dirichlet boundary condition in which the zero modes are not excited and we do not need an IR cutoff. We consider the discrete version on a lattice and use eq. (2.14) to find the corresponding correlators. The results for entanglement entropy are collected in figure 2 where we demonstrate numerical results for different values of dynamical exponent. Based on this figure, it is evident that the value of entanglement entropy increases while the dynamical exponent is increased. Such a behavior is not surprising because as we have seen in the discrete version of the Hamiltonian in eq. (2.5), for larger dynamical exponents the number of lattice points which are coupled together increases and thus the correlation between points inside and outside the entangling region increases. Indeed, for generic z(> 1), the number of correlated points due to the kinetic term increases as 2z + 1 (this is illustrated in any horizontal or vertical slice of different panels of figure 1).

Entanglement measures in Lifshitz-type QFTs
Interestingly, according to figure 2 for large values of z and small enough entangling regions the theory exhibits a volume law scaling of the entanglement entropy. By small enough we mean /N 1. This is in contrast with the behavior of entanglement entropy in typical vacuum states of a local relativistic Hamiltonian which exhibits area law scaling and in (1+1)-dimensions is replaced with a logarithmic scaling. Again, from the expression of lattice Hamiltonian, i.e., eq. (2.5), it is clear that by increasing the critical exponent the theory starts to show nonlocal effects such that for z 1 the corresponding QFT becomes highly nonlocal. This is why we expect and numerically verify the volume law rather than the area law. Actually according to figure 2 one expects that EE is proportional to with a certain power, i.e., S ∼ α . But in order to confirm the volume law, the asymptotic behavior should be α| z→∞ → 1. In the left panel of figure 3 we check this asymptotic behavior considering S S as a function of the dynamical exponent for two different entangling regions. According to these curves, S S approaches to unity while z goes to infinity which confirms the volume law.
Volume law scaling of entanglement entropy has been previously observed in several field theories with nonlocal effects [40][41][42][43][44] (See also [45] for a perturbative study of entanglement entropy in nonlocal theories.). 7 The model of our interest for large dynamical JHEP07(2017)120 exponents is similar to a certain non-local scalar theory studied in [40]. In figure 4 we have shown that this scaling behavior also happens for Renyi entropies for large values of z.
We would like to clarify a very important point about theories which show a volume scaling of entanglement entropy in a certain range of parameters. As long as we are dealing with pure states, one may worry about some fundamental features of entanglement entropy in these theories. Explicitly one may ask about the validity of S A = SĀ (Ā is the complement of A) is such theories for regions which are small enough to show a volume law. Lets consider A the small region which its entropy scales with its volume. In this case the complement regionĀ is a large subregion which dose not have a volume law. If we have the wrong impression that the SĀ scales with its area, then we will lose S A = SĀ! The important point is that SĀ does not scale with the area ofĀ in this case. This is simply because the boundary of A andĀ are the same (∂A = ∂Ā) and the characteristic length of the boundary is smaller than the non-locality scale. In this case it is not enough to consider a strip of dof's around ∂Ā to find the leading term of entanglement entropy. In other words for regions which are large such that their complements show a volume law, the entanglement entropy does not scale with area but with the volume of their complement. With this in mind S A = SĀ is still valid in these theories. For small values of the dynamical exponent we still expect the logarithmic scaling to be the leading contribution to entanglement entropy. For a fixed subregion as the dynamical exponent increases we expect the scaling to tend to volume law. In order to investigate this competition between the contributions due to logarithmic and volume scaling more precisely, we study the ratio of these portions with considering the following fit function 0 is a constant (and also non-universal) term where we drop it in the following discussion by defining a subtracted EE as ∆S (z) ( ) ≡ S (z) ( )−S z 0 . Analysing the numerical data shows that the ratio of the volume term with respect to the subtracted entanglement entropy, i.e., S (z) vol.
∆S (z) should tend to / in the large z limit. This behavior shows that for large values of critical exponents the volume law contribution to EE substantially is dominant. Based on our numerical results, we propose that there should be a leading piece in S (z) ( ), whose dependence on characterizes a transition from area to volume law at a Lifshitz quantum critical point, and its form for d = 1 is given by where for z = 1 should be replaced by a logarithmic scaling. In the right panel of figure 3 we demonstrate the numerical results which are in agreement with our proposal. Of course the scaling introduced in eq. (3.2) is not the most general form but rather the simplest function supporting such a transition between area and volume law. We expect whatever the exact analytic behavior is, which can be found via different methods, e.g., working out the propagator of scalar field on an n-sheeted plane (with Lifshitz symmetry), having the same behavior as eq. (3.2) in the → 0 limit. It is also worth to note that although for large values of z and for N , the scaling of EE differs from the usual area law, but for pure state that we have considered the S A = SĀ condition still holds. We demonstrate this behavior in figure 5.
A relevant question would be how the entropy depends on the dynamical exponent z. We have also studied this numerically which the result is shown in figure 6. As one can see there is a quadratic growth as z starts to increase from the Lorentzian case. The validity of this regime depends on the size of the entangling region. In the large z limit the numerical results shows linear growth of entanglement entropy as a function of z. The smaller the region is, the faster the entropy growth enters the linear regime. This is in agreement with the analytic result found in [36] for a different type of entangling regions called p-alternating sublattice.

Mutual and tripartite information
After analyzing the behavior of entanglement entropy at a critical point with Lifshitz scaling, we now proceed with considering other entanglement measures. In particular we are interested in mutual and tripartite information between different subregions. 9 These quantities are defined in terms of EE as follows where A, B and C are three entangling regions with length a, b and c respectively which we consider in the following configuration: Although the mutual information is always positive, which is a reminiscent of subadditivity of EE, the tripartite information can be positive, negative or zero depending on the quantum field theory and subregion configuration of interest. It is important to mention that in spite of EE which is divergent, these two quantities are defined such that the contributions of different subregions to UV divergent part cancel each other and the resultant expression becomes finite. Using the above definition for mutual and tripartite information and eq. (2.19) one can find these quantities for a Lifshitz QFT.
Let us first focus on mutual information. Here for simplicity we restrict our analysis to a configuration where the length of regions A and B are equal i.e. a = b. The generalization 9 Some aspects of mutual information for z < 1 have been previously studied in [30]. to more generic configurations is straightforward. For such a configuration, we expect the mutual information between two fixed subregions to be a decreasing function of the separation between them which we denote by d. The results are collected in figure 7. Based on this figure there are two points in hand: first, as the dynamical exponent increases, the mutual information between subregions for a fixed configuration increases. This in agreement with what we explained about the increase of EE as the dynamical exponent increases. The second point is related to the well-known fact that mutual information between fixed subregions vanishes while the separations between subregions increases. Our results shows that the critical separation between subregions increases with the dynamical exponent. 10 Once again, this behavior is expected due to the enhancement of spatial correlations for larger values of z. The behavior of tripartite information is qualitatively similar to mutual information. One can check that similar to what happens in Lorentzian free scalar theory, tripartite information is always positive in this Lifshitz scalar theory. This is why we think that our results may not be comparable to some results available in the literature which are found by means of Ryu-Takayanagi proposal on geometries with asymptotic Lifshitz symmetry (see for instance [18,23]). It is worth to remind the reader that in the context of holographic QFTs, it is a well-known fact that holographic tripartite information, at least for states which have a classical gravity dual, must be always negative. The proof of this statement is based on the minimality condition in Ryu-Takayanagi formula [49]. 11

(2 + 1)-dimensions
In this section we proceed our previous analysis in one higher dimension which enables us to more investigate the emergence of volume law and also shape dependence of EE in presence of dynamical exponent z. 12 In order to do so we first generalize our method which was introduced in section 2 to (2 + 1)-dimensions. One can straightforwardly show that JHEP07(2017)120 .
. . . Figure 8. Different entangling regions on a constant time slice of a (2 + 1)-dimensional Lifshitz field theory. From left the first panel shows a disk entangling region and the second one is a square region. The third and forth shapes are deformations of the square region. The area of the shape in the third panel is equal to that of the square but they have different number of corners (6 and 4 respectively). The shape in the forth panel with 8 corners, has the same volume as the square. We call the third one 'a-shape' and the forth one 'v-shape'.
the vacuum correlators with periodic boundary conditions in this case are given by the following expressions where Ω 2 kx,ky = m 2z + 2 sin (N x , N y ) and (k x , k y ) are the number of lattice sites and the momentum in our two spatial directions. It is also straightforward to generalize the vacuum correlators for Dirichlet boundary condition in this case. This is given by Fourier double sine transformation which gives . The rest of the method which we have reviewed briefly in section 2 applies exactly in the same way in (2 + 1)-dimensions. In this section we first work out entanglement entropy for a disk and a square entangling regions (see figure 8). Our main concern here is to verify the emergence of z-dependent scaling behavior of entanglement measures (specifically EE). We also numerically verify that JHEP07(2017)120 The main difference between Lorentzian and Lifshitz cases is the corner contributions to EE are local in the former case but non-local in the latter one. It is well-known in the literature that corner contributions in Lorentzian theories are local effects and thus in (2+1)dimensions which they have a logarithmic behavior, the coefficient of this logarithmic term is an additive function for any number of singularities in the entangling region. Here we present some arguments that the effect of dynamical exponents breaks down the latter property and as z is increased the coefficient of the logarithmic term recedes from being additive.
Disk and square entangling regions. In this subsection we focus on disk and square entangling regions (see the two left panels in figure 8). Once again as we have explained previously (see figure 1), due to the long tail of the kinetic term which causes (2z + 1) lattice sites to be correlated together, we expect that the EE is a monotonically increasing function of z. In figure 9 we have presented the area law behavior of entanglement entropy for square entangling region for z 20. In figure 10 we have shown the behavior of square and disk entangling regions for large values of z. One can see the change in the scaling of EE for regions inscribed in a square with length 20.
We are interested to investigate about the scaling of entanglement entropy for large enough dynamical exponent for different entangling regions. Our numerical data shows that the disk entangling region in (2 + 1)-dimensions for large enough z and small enough regions fits with the following function This function, up its first term which scales with the volume of the disk and is a sign of non-locality in this theory, is the same as what is well-known for disk entanglement in Lorentzian field theories. Our fitting data shows that similar to what we reported in the previous section for (1 + 1)-dimensions, for small regions here also the ratio of the volume term increases while the dynamical exponent increases.
There is one point we would like to clarify here about calculation of entanglement entropy of a disk in (2 + 1)-dimensions on a square lattice. It is sometimes mentioned in JHEP07(2017)120 the literature that entanglement entropy is not a well-defined measure for a disk in (2 + 1)dimensions since it is not a smooth function of the radius of the disk. This is because the EE of a disk shows significant fluctuates if one considers 'continuous' radii for the disks. By continuous we mean radius values which are not multiples of the lattice constant. Due to these fluctuations people focus on mutual information as a more interesting measure for disk entangling region in (2 + 1)-dimensions [1]. 13 Our numerical investigation shows that one can simply get rid of these fluctuations by focusing on disks which their radii is a multiple of the lattice spacing. As one can see there is no such fluctuation in our data in figure 10.
We have also calculated the entanglement entropy of a square entangling region. The main difference between square and disk entangling regions is the singularities (corners) present in the square. It is a well-known fact that these singularities cause new divergent terms in entanglement entropy (see e.g. [53]). Our numerical analysis shows that the entanglement entropy for small squares fits with the following function We postpone a more careful study of the coefficient of this logarithmic term to the next subsection. Here we are interested to investigate how large is the contribution of the volume law to entanglement. Similar to what we reported for the disk, we found that the contribution of the volume law for small regions increases with the dynamical exponent z. Figure 11 is showing the same thing which we previously showed in figure 3 for (1+1)dimensions. Here for a square entangling region with a fixed side , we have plotted the S S in the left panel. This quantity approaches a constant value which supports the volume law behavior in the z → ∞ limit. In the right panel we have plotted the best fit for our proposed function for large values of z.
In figure 12 we have plotted the entropy as function of the volume, i.e. number of sites enclosed by the entangling region. This is done for square and v-shape which have the same volume. One can find that for small regions where the entropy scales with the volume the 13 See also [52] where a F -theorem is introduced in (2 + 1)-dimensions using mutual information for co-centric disks.    Figure 12. Entanglement entropy as a function of number of sites enclosed by the entangling region, n A , for square and v-shape. These two shapes have the same volume (see figure 8). In the left panel we have set z = 20 and in the right panel z = 50. For small regions, where we expect a volume law, the entropy for both shapes coincides. As n A increases, the entanglement of the v-shape overtakes the square since it has a larger area. We have set N x = N y = 250 and m = 0 with Dirichlet boundary condition.
entropy of these two shapes coincides but for larger regions the entropy of the region with a larger area grows faster. For larger values of z the region of coincidence increases. We have also studied the z-dependence of entanglement entropy. In figure 13 we have shown the results for square entangling regions with different size. Our results is similar to what we have found in (1+1)-dimensions. As the dynamical exponent increases from the Lorentzian case, for small values of z the entropy grows quadratically which the length of validity for this regime depends on the size of the entangling region. After this regime the entropy growth enters a linear regime which our results does not show any stop for it.
Our results show that for generic dynamical exponent z, we expect the leading term of entanglement entropy (and also Renyi entropies) to be a function which interplay between these two behaviors. The simplest choice with this property in (2 + 1)-dimensions is Lifshitz theory of interest. In this subsection we will continue our study focusing on the corner contributions to entanglement entropy in (2 + 1)-dimensions. Corner contributions to entanglement entropy are well-known to be local effects. In a local (2 + 1)-dimensional theory corner contributions appear as subleading logarithmic divergent term [1,51,54] S A = S area + S corner log + S 0 . (3.7) where and θ i is the opening angle of the i-th corner. Since the corner contribution is a universal term in the entanglement entropy expansion, a very important question is what kind of characteristic information of the theory it contains. Although we do not know the answer to this question in general, recently in the context of AdS/CFT, 14 it has been shown that in the smooth limit, i.e. θ → π, the coefficient of this term is the same as what appears in the two point function of the stress tensor of the underlying conformal field theory [57]. In order to study these singularity effects in a Lorentzian theory one can consider entangling regions with the same area but different number of singularities. See the two middle panels of figure 8 which show a square (with four corners) and the a-shape which is a deformation of the square which preserves the area but has six corners. In this case we can forget about the leading area term and focus on the subleading logarithmic term which is supposed to capture the difference between these two regions. For these two regions if one works out the entanglement entropy and find the logarithmic term, we expect (from locality) to find which can be verified numerically (see e.g. [53]).
In the case of non-local theories, one would consider two regions with the same volume but different number of corners. In the left panel of figure   square which has the same volume as the square but eight corners (we refer to this shape with index 'v'). Our numerical results for both 'a' and 'v' entangling regions show that as the dynamical exponent increases the ratio of the corner contribution to that of a square starts to decrease which shows that corner contributions are no more local effects. The way we can explain what happens in Lifshitz theory is as follows: contribution of each corner is mixed up with others due to correlation of (2z + 1) lattice cites. This has been demonstrated in figure 14. In order to have a neat picture we have just focused on the nearest points to the singularity which are show in green. This is done for 'v'-shape and the similar thing happens for the 'a'-shape region. The red shadowed points are those which are in the (2z + 1) neighbourhood of the green points. The left panel belongs to the Lorentzian case z = 1 which non of the shadowed points of each green point coincide with each other. The right panel belongs to z = 6. In this case there are lots of shadowed points which belong to the shadow of more than one corner (green) point. The number of such points increases with the dynamical exponent and this is why we think the corner contributions are no more additive in the case of Lifshitz scalar theory. This mixing procedure will also happen in Lorentzian theories if there are adjacent corners in the entangling region [58].

Entanglement measures in massive scalar theory
In this section we continue our study by considering the massive scalar field theory. In order to investigate the effects of nonzero mass on the entanglement we will consider the periodic boundary condition for the field in different dimensions. The results for the entanglement and Renyi entropy in (1 + 1)-dimensions are summarized in figure 15.
According to these results similar to the massless case the EE increases with the dynamical exponent. The important thing to note here is that the massive case is completely different from the massless case. We have not found volume law regime for small subregions at any value of z in our study. Our results shows that for periodic boundary conditions, for mass parameters which lead to stable numerical results, we always recover an area law scaling for entanglement entropy. Increasing the mass parameter makes the entanglement to grow slower as a function of the subregion length and also the dynamical exponent. We have plotted some results for entanglement and Renyi entropies in figure 15. In case of Dirichlet boundary condition for small mass parameters in the m → 0 limit we recover the results reported in the massless section. We have not found any range of parameters even with large dynamical exponents with a volume law behavior for small regions. This is in agreement with the behavior of the correlator in (1 + 1)-dimensions in the small r limit which is given by for integer dynamical exponent z > 3. The scaling of the leading term does not depend on the dynamical exponent and is divergent in the m → 0 limit. This causes our intuitive picture (illustrated in figure 1) not to be valid in the massive regime. A more careful analysis is needed to understand the massive regime of this theory which we postpone it for future works.

Conclusions and discussions
We have mainly studied the entanglement structure of a free massless scalar theory (see equation (2.1)) with Lifshitz scaling symmetry. This family of theories are supposed to describe systems at quantum critical points. We have analysed our theory of interest in 1+1 and 2+1 space-time dimensions. We explored the role of the dynamical exponent in entanglement measures. We would like to first summarize our main results: • In the massless limit for small entangling regions the leading part in the EE expansion transits from area law to volume law for large dynamical exponents. Regarding to our results in 1 and 2 spatial dimensions for several entangling regions, we propose the following scaling behavior as the leading term of EE in a (d + 1)-dimensional theory where is a characteristic length of the entangling region. Indeed, the spatial part of the kinetic term in this theory couples (2z + 1) lattice sites together at each point, and therefore the value of entanglement grows due to the correlation between more JHEP07(2017)120 points in comparison with the Lorentzian case. In other words the theory shows more non-local effects for larger dynamical exponents. This is why one would naturally expect a volume law scaling rather than an area law for small entangling regions.
Regarding to this result we would like to also note that this result is in contrast with what has been previously reported in [16] where the author has argued that in any theory with a generic form of dispersion relation as ω = f (k 2 ), where k µ is the momentum vector, the EE scales with the area of the entangling region.
• We have shown that the mutual information is an increasing function of the dynamical exponent z. Also the fall off mutual information in large separation limit is slower in comparing to the relativistic case. This behavior is also due to stronger correlations between lattice points for larger values of z.
• Similar to free scalar Lorentzian theory, the resultant tripartite information is always positive for all values of z. This shows that at least the vacuum state of this theory can not have a solution of a classical gravity theory as a holographic dual. This is the main clue that forbids us to compare our results (especially eq. (4.1)) with some previous studies on Lifshitz theories in the context of gauge/gravity duality.
• In 2 spatial dimensions we have argued that corner contributes which are known to be local effects in Lorentzian field theories become more and more non-local while the dynamical exponent increases. As a result corner contribution to entanglement entropy are no more additive in these theories.
We have also observed that for massive theory the behavior of the m → 0 limit at fixed z is very similar to the behavior of large z at fixed m. Although at the moment it is not an easy task to understand this behavior from the correlation functions for generic z because of technical reasons, but it is worth to note that this behavior is in agreement with what is known about a massive scalar field propagating on a Lifshitz background in the context of holography (see eq. 3.23 of [24]).
Of course it would be of great interest if one can verify our results with a concrete analytical analysis. As an example one can study propagators of free scalar theories on a n-sheeted space-time and work out the leading behavior entanglement and Renyi entropies analytically. Recently a nice paper appeared where heat kernels for non-relativistic field theories on specific manifolds has been worked out in it [60,61]. In principle it is possible to extend this study to singular manifolds and work out the entanglement entropy.
Finally we would like to note some interesting questions which generate new directions for our analysis. One way to extend our study is to go beyond the vacuum state. In the case of mixed state, such as thermal states it is interesting to analyse other entanglement measures which are more suitable for mixed state, e.g., logarithmic negativity [7,8,59].
Another interesting question is about the time evolution of EE in such theories which is now a typical way to study universal behaviors of these family of theories in out of equilibrium phases. One can study entanglement measures under quantum quenches and study the corresponding thermalization process. Indeed, this analysis may shed light on JHEP07(2017)120 the role of critical exponent on the entanglement structure. We will report some interesting results in these directions in early future [62].