Entanglement gap, corners, and symmetry breaking

We investigate the finite-size scaling of the lowest entanglement gap $\delta\xi$ in the ordered phase of the two-dimensional quantum spherical model (QSM). The entanglement gap decays as $\delta\xi=\Omega/\sqrt{L\ln(L)}$. This in contrast with the purely logarithmic behaviour as $\delta\xi=\pi^2/\ln(L)$ at the critical point. The faster decay in the ordered phase reflects the presence of magnetic order. We analytically determine the constant $\Omega$, which depends on the low-energy part of the model dispersion and on the geometry of the bipartition. In particular, we are able to compute the corner contribution to $\Omega$, at least for the case of a square corner.


Introduction
In recent years the cross-fertilization between condensed matter and quantum information fueled an impressive progress in our understanding of quantum many-body systems [1,2,3,4]. The entanglement spectrum (ES) has been the subject of intense investigation. Let us consider a system in its ground state Ψ⟩ and a spatial bipartition of it as A ∪Ā (see Fig. 1). The reduced density matrix ρ A = TrĀ Ψ⟩⟨Ψ of A can be written as Here H A is the so-called entanglement hamiltonian. The entanglement spectrum levels ξ i = − ln(λ i ), with λ i the eigenvalues of ρ A , are the "energies" of H A . Early works [5,6,7] on entanglement spectra aimed at understanding the effectiveness of the density matrix renormalisation group (DMRG) [8,9] to simulate one-dimensional systems.
Recently, an intense theoretical activity has been devoted to understand the ES in fractional quantum Hall systems [10,11,12,13,14,15,16,17,18,19,20,21], topologically ordered systems [22,23,24], magnetically ordered systems [25,26,27,28,29,30,31,32,33,34,21,35,36,37], Conformal Field Theories (CFTs) [38,39,40,41], and systems with impurities [42]. The entanglement gap (or Schmidt gap) δξ emerged as a natural quantity to investigate. δξ is the gap of the entanglement hamiltonian, and it is defined as where ξ 0 and ξ 1 are the first two low-laying ES levels. For the standard energy gap, i.e., the gap of the physical hamiltonian, there exists a "universal" correspondence between its scaling behaviour and ground state properties, such as the decay of correlation functions [43]. Much less is known for the entanglement gap, although several results are available. For instance, its behaviour at one-dimensional quantum critical points has been investigated [44,45,5,6,28,11,27,32,46]. In CFTs it is well established that δξ decays logarithmically as δξ ∝ 1 ln( ) with the subsystem length [38]. Similar scaling is found in models that are solvable via the corner transfer matrix technique [44]. Higher-dimensional models are uncharted territory. Interestingly, some explicit counterexamples show that the closure of the entanglement gap does not signal criticality [21], also for the momentum-space ES [47]. The scenario is different deep in ordered phases of matter. For instance, the lower part of the ES of magnetically-ordered ground states that break a continuous symmetry [29] is reminiscent of the Anderson tower-of-states [48,49,50]. This has been verified in systems of quantum rotors [29], in the two-dimensional Bose-Hubbard model in the superfluid phase [31] (see also [37]), and also in Heisenberg antiferromagnets on the square [34] and on the kagome lattice [36].
In the tower-of-states scenario gaps in the lower part of the ES decay as a power-law with the subsystem volume, with multiplicative logarithmic corrections [29]. Higher ES gaps exhibit a slower decay [29,31,35].
Given the lack of general results, exactly solvable models can provide valuable insights into the generic features of the entanglement gap. Here we investigate the entanglement gap in the ordered phase of the two-dimensional quantum spherical model (QSM) [51,52,53,54,55]. Despite its appealing simplicity, the QSM contains several salient features of generic quantum many-body systems. The model is mappable to a system of free bosons with an external constraint, implying that its properties can be studied with moderate cost. Its classical version proved to be valuable to validate the theory of critical phenomena and finite size scaling [56]. The ground state phase diagram of the two-dimensional QSM exhibits a paramagnetic (disordered) phase and a ferromagnetic (ordered) one, which are divided by a continuous quantum phase transition. The universality class is that of the three-dimensional classical O(N ) vector model [57] in the large N limit [58,52,53]. Entanglement properties of O(N ) models have been addressed in the past [59,60] (see also [61,62,63,64] for recent studies in the QSM). Here we consider a two-dimensional lattice of linear size L. The typical bipartitions that we use are reported in Fig. 1. Figure 1 (a) shows a bipartition with a straight boundary. The bipartition in Fig. 1 (b) contains a square corner. The effects of corners in the scaling of entanglement-related quantities is nontrivial, and it has been studied intensely in the last decade [65,66,67,68,69,70,71,4,72].
Our main result is that in the ordered phase of the QSM, in the limit L, x , y → ∞ with the ratios ω x,y = x,y L (see Fig. 1) fixed the entanglement gap decays as Here the dots denote subleading terms that we neglect. The constant Ω, which we determine analytically, depends on the low-energy properties of the model and on the geometry of the bipartition. In particular, we analytically determine the corner contribution to Ω. The "fast", i.e., power-law behaviour as 1 √ L in (3) reflects the presence of magnetic order, whereas the logarithmic correction is similar to the critical behaviour [64] of δξ. Finally, we should mention that Eq. (3) is different from the result derived in Ref. [29], where it was shown that for O(N ) models δξ ∝ (L ln(L)) −1 .
The manuscript is organised as follows. In section 2 we introduce the QSM. In section 3 we review the finite-size scaling of the ground-state two-point correlation functions. In section 4 we briefly overview the calculation of the entanglement gap in the QSM. Section 5 is devoted to the derivation of our main results. In section 6 we provide numerical checks. We conclude in section 7. In Appendix A we report some technical derivations.

Quantum Spherical Model
The quantum spherical model [52,53,54] (QSM) on a two dimensional cubic lattice of volume V = L 2 , with L being the lattice linear size, is defined the hamiltonian Here, n = (n x , n y ) denotes a generic lattice site, and ⟨n, m⟩ a lattice bond joining two nearest-neighbour sites. We set J = 1 in (4). The spin s i and momenta p i variables satisfy standard bosonic commutation relations Here we refer to the parameter g as the quantum coupling. Indeed, in the limit g → 0 the model reduces to the classical spherical model [73,74]. The spherical parameter µ is a Lagrange multiplier that fixes the global magnetization as n ⟨s 2 n ⟩ = V, To diagonalize the QSM hamiltonian (4), one can exploit its translational invariance. First, one performs a Fourier transform as where the sum is over k = (k x , k y ) in the first Brillouin zone k i = 2π Lj, with j ∈ [−L 2, L 2] integer. In Fourier space one obtains The single-particle dispersion relation is given as To fully diagonalise (8) we introduce the new bosonic ladder operators b k and b † k as where α 2 k = g 2Λ −1 k . Now, the hamiltonian (8) is fully diagonal, and it reads as For the following, it is useful to consider the ground-state two-point correlation functions ⟨s n s m ⟩ and ⟨p n p m ⟩. They are given as [54] S nm = ⟨s n s m ⟩ = Importantly, the trivial identity holds By using (12), the spherical constraint (6) can be rewritten as Eq. (16) is the so-called gap equation in the context of the large-N model [75]. A crucial observation is that the correlator (12) exhibits a singularity for k = 0. This zero mode will play a crucial role in the behaviour of the entanglement gap.
In two dimensions at zero temperature the QSM exhibits a second-order phase transition at a critical value g c . The value of g c is known analytically as For g < g c the QSM exhibits a magnetically ordered phase, which is the focus of this work. At g > g c the ground state is paramagnetic. Different phases are associated with different behaviour of the spherical parameter µ. In the paramagnetic phase one has that µ is finite nonzero. On the other hand, µ = 0 at critical point, and in the ordered phase. The different phases of the model correspond to different finite-size scaling behaviours of µ. In the paramagnetic phase one has µ = O(1). At the critical point one can show that µ = O(1 L 2 ). In the ordered phase µ = O(1 L 4 ). The critical behaviour at g c is in the universality class of the large-N N -vector model in three dimensions [53].

Spin and momentum correlators
Here we summarise the finite-size scaling of the spin-spin correlation function S nm (cf. (12)) and the momentum correlation function P nm (cf. (13)). Let us focus first on the spin correlator. We are interested in the limit L → ∞. We can decompose the correlator as The first term is the leading term in the large L limit. Note that the first term depends on L because µ depends on L. The second term in (18) is the first subleading in powers of 1 L. The dots denote more subleading terms that we neglect. The thermodynamic contribution is given as The finite-size part has the surprisingly simple form [64] as Here we defined The prime in the sum means that one has to remove the term with (l, l ′ ) = (0, 0). Eq. (20) holds in the limit L → ∞ and µ → 0, i.e., for g ≤ g c . The correlators S nm depend only on n x − m x and n y − m y , reflecting translation invariance. Moreover, the infinite sums in l, l ′ enforces that S nm is periodic along the two directions, i.e., it is invariant under n y − m y → n y − m y ± L and n Fig. 1 (a)), whereas no singularity occurs for ω x < 1 and ω y < 1, i.e, in the presence of a bipartition with a corner (see Fig. 1 (b)). Let us consider the case ω y = 1. Now the terms with l = 0 and l ′ = ±1 in (20) are singular in the limit n x − m x → 0 and n y − m y → ±1. Terms with l ′ > 1 or l > 1 in (20) are not singular, and do not affect the singularity structure of S nm . These singularity will play an important role in section 5. Similar to (18), we can decompose the momentum correlator as Here we defined The finite-size part P nm has the same structure as (20), and it reads as with F ll ′ (n, m) as defined in (21). Eq. (23) and Eq. (24) are obtained from Eq. (19) and Eq. (20) by using (15). As for Eq. (20), the finite-size term (24) is singular if subsystem A spans the full lattice in one of the two directions, i.e., if ω y = 1 or ω y = 1 (see Fig. 1).
For ω y = 1 the singularity occurs for l = 0 and l ′ = ±1 in the limit n x − m x → 0 and n y − m y → ±1. Finally, the first term has a stronger singularity than the second one.

Spherical parameter
Let us discuss the finite-size scaling of the spherical constraint µ (cf. (16)) in the ordered phase of the QSM. For g ≤ g c the spherical parameter vanishes in the thermodynamic limit. At the critical point one has the behaviour µ ∝ γ 2 2 (2L 2 ), with γ 2 a universal constant. To derive the behaviour of µ in the ordered phase we use Eq. (12) in the gap equation (16). We obtain As it clear from the exponent in the last term in (25) the scaling as µ ∝ 1 L 2 at criticality implies that terms with large l, l ′ are exponentially suppressed. On the other hand, for µ ∝ 1 L 4 this is not the case because the term in the exponent in (25) is O(1 L). First, we anticipate that the second term in (25) is O(1 L 2 ), and it is subleading. To extract the leading behaviour of µ we use the elegant identity involving the function K(σ, d, y) defined as [76] K(σ, d, y) = ′ l(d) Here the sum is over the d-dimensional vector of integers l i ∈ (−∞, ∞), K σ (z) is the modified Bessel function of the second kind, and y > 0 and σ are real parameters. We are interested in the case d = 2 and σ = 1 2 (cf. (25)). One can show that [76] The constant C(σ, d) for d = 2 reads as where ζ(x) is the Riemann zeta function, and β(x) is the analytic continuation of the Dirichlet series To apply (27) we fix y = µ 2L. In the limit µ → 0 the leading behaviour of K is given by the first term on the right hand side in (27). After using that in (25) we obtain In (30) we are neglecting vanishing terms in the limit L → ∞. The second term in (30) is also simply obtained by isolating the term with k = 0, i.e., the zero mode, in the sum in (16). It is now clear that we can parametrize µ as After substituting in (30), we obtain that Note that the constant γ 4 is not universal, as it is clear from the explicit dependence on g. This is expected, and it is in contrast with the result at the critical point, where µ = γ 2 2 (2L 2 ), with γ 2 universal.

Entanglement gap in the QSM
Here we briefly review how to calculate the entanglement gap in the QSM. Entanglement properties of the QSM are derived from the two-point correlation functions (12) and (13) because the model can be mapped to free bosons (see Ref. [7] for a review). We first define the correlation matrix C restricted to subsystem A as with S A and P A defined in (12) and (13), with n, m ∈ A. Since the QSM is mapped to free-bosons, the reduced density matrix of a subsystem A is a quadratic operator, and it is written as [7] Here H A is the so-called entanglement hamiltonian, k are single-particle entanglement spectrum levels, and b k are free-bosonic operators. Z is a normalization factor. The eigenvalues e k of C A are obtained from the k as The entanglement spectrum, i.e., the spectrum of H A is obtained by filling in all the possible ways the single-particle levels k . The lowest ES level is the vacuum state. Thus, the lowest entanglement gap δξ (Schmidt gap) is with 1 the smallest single-particle ES level, or equivalently, the largest e 1 (cf. Eq. (35)).

Scaling of the entanglement gap in the ordered phase of the QSM
In this section we investigate the scaling of the entanglement gap for g < g c , i.e., in the ordered phase of the QSM. First, it has been numerically observed in Ref. [64] that for g < g c , in the limit L → ∞ the flat vector 1⟩ defined as with A = x y , is the right eigenvector of C A corresponding to the largest eigenvalue e 1 , i.e., the zero-mode eigenvector. Moreover, 1⟩ is also eigenvector of the matrix S A . It is interesting to investigate the structure of the associated eigenvalue. This is calculated as After using (18), it is straightforward to numerically check that the thermodynamic part of the correlator S (th) for large L gives a subleading term as L ln(L) in (38) (see section 6). The leading contribution is given by the finite-size part of the correlator S (L) , and it is O(L 2 ). An important observation is that due to the scaling as µ = γ 2 4 L 4 , the dependence on the coordinates n, m in (20) can be neglected. Thus, a straightforward calculation yields One should observe that Eq. (39) is exactly the contribution of k = 0 in the sum in (12).
Physically, this means that in the ordered phase of the QSM for g < g c the leading behaviour of the eigenvalue of S A associated with the zero-mode is simply obtained by isolating the term with k = 0 in (12). This happens because of the "fast" decay as µ ∝ 1 L 4 . This is not the case at the critical point [64], where µ ∝ 1 L 2 . Moreover, this result suggests that one can decompose the correlator S as Here s 0 ∝ L 2 , and the dots are subleading terms that we neglect. By using (40) and the fact that P is finite in the limit L → ∞, it is straightforward to show that the eigenvalue e 1 of C A = P ⋅ S in the limit L → ∞ is given as (see [77] and [64]) Here we have To proceed we now use that the expectation value ⟨1 P 1⟩ decays as ln(L) L, i.e., with a multiplicative logarithmic correction. Note that the same scaling behaviour is observed at the critical point [64]. The derivation requires minimal modifications as compared with the critical case, and it is reported in Appendix A. The main ingredients are standard tools in the finite-size scaling theory, such as Poisson's summation formula and the Euler-Maclaurin formula. Let us discuss the final result. Clearly, we can treat the contribution of the thermodynamic part (cf. (23)) and the finite-size part (cf. (24)) separately. Similar to what happens at the critical point [64], the finite-size part contributes only if the boundary between the two subsystems is straight. For simplicity we consider the bipartition with ω x = 1 p and ω y = 1 q, with p, q ∈ N. For ω y < 1 the boundary between the two subsystems is not straight, i.e., it has square corner. One obtains the thermodynamic contribution as The function η p ′ ,q ′ (k x , k y ) reads as
The dots in the brackets denote terms with higher powers of 1 (k x + p ′ p). These can be derived systematically by using the Euler-Maclaurin formula. The function ψ ′ (x) is the first derivative of the digamma function ψ(x) with respect to x. The behaviour as ln(L) L is clearly visible in (44). Similar to the critical point [64], η p ′ ,q ′ is determined by the low-energy part of the dispersion of the QSM. Finally, let us consider the finitesize contribution (24). From (24) it is clear that the finite-size correlator is regular for ω y < 1 and ω x < 1, whereas it exhibits a singularity for ω y = 1 or ω x = 1, i.e., for the case of straight boundary (see Fig 1 (b)). Again, this is exactly the same at the critical point [64]. The derivation of the singular contribution is reported in Appendix A.2. The final result reads The minus sign in (45) where the constant Ω ′ encodes information about the geometry of the bipartition and the model dispersion. In Eq. (46) we neglect subleading terms in the limit L → ∞. From (46), after using (35) one obtains that Here Ω ′ is defined in (46), and we neglect subleading terms. Few comments are in order. First, in the ordered phase δξ vanishes in the thermodynamic limit as a power law with L, except for a logarithmic correction. This is different at the critical point, where [64] 100 1000 g=5 ω y =1 ω x =1/2 (numerical) analytical Figure 4. Rescaled flat-vector expctation value ⟨1 P 1⟩L of the momentum operator in the ordered phase of the QSM. We show data for several bipartitions with aspect ratios ω x , ω y (see Fig. 1). For ω y < 1 the boundary between the two subsystems is not smooth (see Fig. 1 (b)). Symbols are exact numerical results. The dashed-dotted lines are analytic predictions from (43) and (45). δξ = π 2 ln(L). The power-law decay of the entanglement gap in symmetry-broken phases has been also numerically observed in magnetic spin systems [34,36] and in the ordered phase of the two-dimensional Bose Hubbard model [31]. Note, however, that even with state-of-the-art numerical methods it is challenging to observe the logarithmic correction. Finally, in Ref. [29] it has been suggested that in the presence of continuous symmetry breaking the gaps in the lower part of the entanglement spectrum are δξ ∝ (L ln(L)) −1 .
This different from (47) (note the square root in (47)). The unexpected square root in Eq. (47) could be explained by the way in which in the QSM the spherical constraint is enforced (cf. (16)). Further study would be needed to clarify this issue. Finally, it is interesting to understand the behaviour of δξ as the critical point is approached from the ordered side of the transition. A natural scenario is that upon approaching the transition the 1 √ L is "gapped" out and it gives an extra 1 ln(L), which allows to recover the expected result [64] δξ ∝ 1 ln(L).

Numerical results
In this section we provide numerical evidence supporting the analytic result derived in section 5. Let us start discussing the finite-size scaling of the expectation value ⟨1 S 1⟩. We report numerical data in Fig. 5, for fixed g = 3 (circles) and g = 5 (squares). We only show data for the bipartition with straight boundary ω y = 1 (see Figure 1  g=5 ω y =1 ω x =1/2 (numerical) g=5 ω y =1/2 ω x =1/4 (numerical) analytical Figure 5. Largest eigenvalue e 1 of the correlation matrix C restricted to A. We plot e 1 L versus L. Note the logarithmic scale on the x-axis. Symbols are exact numerical data. The dashed-dotted lines are analytic predictions. Note that for ω y the boundary between the two subsystems has a corner (see Fig. 1 (b)).
in the figure is the analytic result in Eq. (39), which is in perfect agreement with the numerical data. Again, we should stress that Eq. (39) originates only from the finitesize part S (L) (cf. (18)). Thus, it is interesting to investigate the finite-size scaling of the flat-vector expectation value calculated using the thermodynamic contribution S (th) . We report this analysis in the inset of Fig. 5 plotting ⟨1 S (th) 1⟩ L versus L. Data are for g = 5. Interestingly, the figure shows that ⟨1 S (th) 1⟩ ∝ L ln(L). This confirms that at the leading order in L the expectation value ⟨1 S 1⟩ is dominated by the contribution of the zero mode. Finally, we should mention that it would be interesting to clarify the origin of the logarithmic divergence of the thermodynamic contribution. Let us now discuss the flat-vector expectation value of the momentum correlator ⟨1 P 1⟩. In contrast with the spin correlator, both the thermodynamic and the finite-size part (cf. (22)) contribute to the leading behaviour at large L. Our numerical data are reported in Fig. 6. In the figure we plot ⟨1 P 1⟩L versus L. We show data for ω x = 1 2, ω y = 1 and ω y = 1 2. Note that for ω y = 1 the boundary between A and its complement is straight. The numerical data in Fig. 6 confirm the expected behaviour as ln(L) L in Eq. (43) and Eq. (45). For ω y = 1 the prefactor of the logarithm is obtained by summing Eq. (43) and Eq. (45), whereas in the presence of a square corner only Eq. (43) has to be considered. Finally, we discuss the largest eigenvalue e 1 of the restricted correlation matrix C A (cf. (33)). The entanglement gap δξ is obtained from e 1 via Eq. (35). Note that the vanishing of δξ is reflected in a diverging e 1 . We show numerical data for e 1 L in Fig. 6 plotted versus L. We consider several aspect ratios ω x and ω y , focusing on g = 5. In all the cases the data exhibit the behaviour e 1 = Ω ′ L ln(L). The constant Ω ′ , which depends on the geometry and on low-energy properties of the QSM is obtained by combining Eq. (41) with Eq. (39) (43) (45). The analytic predictions are reported in Fig. 6 as dashed-dotted lines and are in perfect agreement with the numerical data. This implies that the entanglement gap δξ satisfies (47).

Conclusions
We investigated the entanglement gap in the magnetically ordered phase of the twodimensional QSM. Our main result is that the entanglement gap decays as δξ = Ω L ln(L). We analytically determined the constant Ω, which depends on the geometry of the bipartition and on the low-energy physics of the model.
There several intriguing directions for future work. First, it would be interesting to explore whether is possible to extend our results to the N -vector model at finite N . An interesting question is whether the discrepancy with the results of Ref. [29] can be attributed to the large N limit. Furthermore, it is important to understand how the scaling of the entanglement gap depends on dimensionality. This issue could be easily addressed because the QSM is exactly solvable in any dimension. Another intriguing direction is to further investigate the role of corners. For instance, it would be interesting to investigate the dependence of the entanglement gap on the corner angle. It would be also interesting to investigate how the outlined scenario is affected by long-range interactions. This should be straightforward because the QSM is exactly solvable also in the presence of long-range interactions. An exciting possibility is to investigate what happens to the entanglement gap in the presence of disorder [78,79,80,81]. Finally, a very interesting direction is to study δξ after a quantum quench. This could be addressed using the results of Ref. [82].

Acknowledgements
I would like to thank Sascha Wald and Raul Arias for several discussions in a related project. I acknowledge support from the European Research Council under ERC Advanced grant 743032 DYNAMINT.
Appendix A. Derivation of the expectation value ⟨1 P 1⟩ In this appendix we derive the large L behaviour of the expectation value of the momentum correlator with the flat vector ⟨1 P 1⟩ (cf. (42)). We consider the leading, i.e, the thermodynamic limit, as well as the first subleading contributions. The main goal is to show that the expectation value exhibits a multiplicative logarithmic corrections. Two types of contributions are present. One originates from the thermodynamic limit of the correlator, whereas the second one is due to the first subleading. The latter is present only for straight boundary between the two subsystems (see Fig. 1).

Appendix A.1. Thermodynamic contribution
Here derive the thermodynamic contribution, which is given as ⟨1 P (th) 1⟩. Here 1⟩ is the flat vector in region A, i.e, The momentum correlation reads After performing the sum over n and m in (A.2), and after changing variables to k ′ x = Lω x k x π and k ′ y = Lω y k y π, we obtain To extract the large L behaviour of (A.3) it is useful to split the integration domains [0, Lω x 2] and [0, Lω y 2] to write We now restrict ourselves to the case with ω x = 1 p and ω y = 1 q, with p, q positive integers. After a simple shift of the integration variables as k x → k x − l x ω x and k y → k y − l y ω y , one obtains We now focus on the behaviour at the quantum phase transition. We set g = g c , µ = γ 4 L 4 (cf. (31)), and we expand (A.5) in the limit L → ∞. This gives Here γ 4 is defined in (31), and we used the periodicity of the trigonometric functions. The term γ 4 L 2 can be neglected for L → ∞. Importantly, as a result of the large L limit, Eq. (A.6) depends only on the low-energy part of the dispersion of the QSM, although it contains non-universal information. We now have to determine the asymptotic behaviour of the sum over l x , l y in (A.6), i.e., of the function η p ′ ,q ′ (k x , k y ) defined as The asymptotic behaviour of η in the limit L → ∞ can be obtained by using the Euler-Maclaurin formula. Given a function f (x) this is stated as In (A.8) the dots denote terms with higher derivatives of f (x) calculated at the integration boundaries x 1 and x 2 . These can be derived to arbitrary order. To proceed, we first isolate the term with either l x = 0 or l y = 0 in (A.7). The remaining sum after fixing l x = 0 or l y = 0 can be treated with (A.8). We define this contribution to the large L behaviour of η p ′ ,q ′ as η 0 , which is given as In the derivation of (A.9) we neglected the boundary terms in (A.8) because they are subleading.
We are now left with the sums over l x ∈ [1, L (2p)] and l y ∈ [1, L (2q)] in (A.7). These be calculated again by using (A.8). We first apply (A.8) to the sum over l x . We have two contributions. The first one is obtained after evaluating the integral in (A.8) at x 2 = L (2p). After expanding the result for L → ∞, we obtain the contribution η 1 given as Note the term ln(L) L in (A.10). The sum over l y in (A.10) can be performed exactly to obtain in the large L limit, yielding Here ψ ′ (z) is the first derivative of the digamma function ψ(z) with respect to its argument. The second contribution is obtained by evaluating the integral in (A.8) at x 1 = 1. The remaining sum over l y cannot be evaluated analytically. However, one can, again, treat the sum over l y with (A.8). After neglecting the boundary terms in (A.8), which are subleading for large L, and after evaluating the integral in (A.8) at x 2 = L (2q), we obtain the contribution η 2 as To obtain the full contribution of the sum over l x in (A.7) we now have to consider the effect of the boundary terms in (A.8).
Having discussed the contribution which derives from approximating the sum over l x in (A.7) with the integral in (A.8), we finally focus on effect of the boundary terms in (A.8). Let us consider the first boundary term (first term in the second row in (A.8)). A term as ln(L) L is obtained by fixing l x = 1, other contributions are subleading. After performing the sum over l y one obtains the first boundary contribution η b1 as In a similar way the second boundary term (last term in (A.8)) gives (A.14) Note that boundary terms in (A.8) are expected to be small. Specifically, the k-th term is suppressed as 1 (k + 1)!. The final result for η(k x , k y , p, p ′ , q, q ′ ) is obtained by putting together (A.9)(A.11) (A.12)(A.13)(A.14) to obtain η p ′ ,q ′ (k x , k y ) = η 0 + η 1 + η 2 + η b1 + η b2 . In this section we derive the leading behaviour in the large L limit of ⟨1 P (L) 1⟩.
Interestingly, we show that in the presence of a straight boundary between the two subsystems (see Fig. 1) one has the behaviour ⟨1 P (L) 1⟩ ∝ ln(L) L. On the other hand, in the presence of corners the multiplicative logarithmic correction is absent. The finite-size correlator reads as (cf. (24) Crucially, if ω x < 1 and ω y < 1, the denominators in (A.16) are never singular. This implies that the logarithmic correction is not present, which can be straightforwardly checked numerically. Let us now consider the situation with ω x < 1 and ω y = 1. Now, a singularity appears in the limit L → ∞ for l = 0 and l ′ = ±1. We numerically observe that only the first term in (A.16) gives rise to a singular behaviour. Thus, we neglect the second term and fix l = 0, obtaining To proceed, let us now consider the case with l = 1. It is clear that the contribution from l = −1 is the same. We can restrict the sum over x in (A. 19) to x > 0 because of the symmetry x → −x. We also restrict to y < 0 because the singularity in (A.19) occurs for y < 0. We now have where the contribution of l = −1 in (A.19) has been taken into account.