Free energy asymptotics of the quantum Heisenberg spin chain

We consider the ferromagnetic quantum Heisenberg model in one dimension, for any spin \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$S\ge 1/2$$\end{document}S≥1/2. We give upper and lower bounds on the free energy, proving that at low temperature it is asymptotically equal to the one of an ideal Bose gas of magnons, as predicted by the spin-wave approximation. The trial state used in the upper bound yields an analogous estimate also in the case of two spatial dimensions, which is believed to be sharp at low temperature.


Introduction
The ferromagnetic quantum Heisenberg model is one of the most important and widely studied models of statistical mechanics. In dimensions d ≥ 3, the model is widely believed to display long-range order at low temperature, but a rigorous proof remains elusive. Based on the concept of long-range order, the low temperature properties of the model are usually examined using spin-wave theory. In the spin-wave approximation, one assumes that the low-energy behavior of the system can be described in terms of collective excitations of spins called spin waves. From an equivalent point of view, which dates back to Holstein and Primakoff [17], these spin waves are known as bosonic quasiparticles called magnons.
The spin-wave approximation has been very successful, predicting for example a phase transition in three and more dimensions, or the T 3/2 Bloch magnetization law [7,8]. In his seminal 1956 paper [14], Dyson derived further properties of the quantum Heisenberg model which, among other things, included the low temperature expansion of the magnetization.
While there was little doubt about the validity of spin-wave theory in three (or more) dimensions, a rigorous proof of some of its predictions has only recently been given in [13] (see also [12]). There it was proved that the free energy of the three-dimensional ferromagnetic quantum Heisenberg model is, to leading order, indeed given by the expression derived using spin-wave approximation, for any spin S ≥ 1/2 (see also [10,25] for earlier non-sharp upper bounds, or [5,11] for results in the large S limit).
The situation is different in lower dimensions. It has been known since the seminal work of Mermin and Wagner [19] that the d = 1 and d = 2 dimensional quantum Heisenberg models do not exhibit long-range order at any nonzero temperature. The low temperature behavior of the system in low dimensions is thus very different from the one in three or higher dimensions, and it is less clear whether spin-wave theory should also be valid in lower dimensions.
In 1971, Takahashi [22] derived a free energy expansion for d = 1 in the case S = 1/2. In this special case, the quantum Heisenberg model is exactly solvable via the Bethe ansatz [6]. The spectrum of the (finite size) model can be obtained by solving the corresponding Bethe equations. Under certain assumptions (known as string hypothesis) on the solutions of these equations, he derived what are now known as thermodynamic Bethe equations, an analysis of which leads to a formula for the free energy. Later, in [23] he derived an alternative free energy expansion using (a modified) spin-wave theory (for any S, and also in two dimensions). Interestingly, the second terms in the (low temperature) free energy expansions in [22,23] do not agree with the predictions of conventional spin-wave theory [7,8,14,17]. (The leading terms do agree, however.) The thermodynamic Bethe equations have been used not only for the Heisenberg spin chain, but also in other models including the Kondo model [1][2][3]21] or the Gross-Neveu model in high energy physics [4]. For more applications of the string hypothesis and its relation to numerous other models in physics, we refer to the review articles [18,24].
In the present paper, using different methods, we prove that, to leading order, the formula derived by Takahashi based on the Bethe ansatz and the string hypothesis in [22] is indeed correct. Our analysis does not use the Bethe ansatz and our result holds for any spin S. It therefore also partly justifies the spin-wave approximation derived in [23]. We shall utilize some of the methods developed for the three-dimensional case in [13], but novel ingredients are needed to treat the case of lower dimensions, both for the upper and the lower bounds.

Model and main result
We consider the one-dimensional ferromagnetic quantum Heisenberg model with nearest neighbor interactions. For a chain of length L, it is defined in terms of the Hamiltonian (2.1) Here, S = (S 1 , S 2 , S 3 ) denote the three components of the spin operators corresponding to spin S, i.e., they are the generators of the rotations in a 2S + 1-dimensional representation of SU (2). The Hamiltonian H L acts on the Hilbert space H L = L x=1 C 2S+1 . We added a constant S 2 for every bond in order to normalize the ground state energy of H L to zero.
Our main object of study is the specific free energy for β > 0, and its thermodynamic limit We are interested in the behavior of f (S, β) in the low temperature limit β → ∞ for fixed S. Our main result is as follows.
where ζ denotes the Riemann zeta function.
The proof of Theorem 2.1 will be given in Sects. 4 and 5, where we derive quantitative upper and lower bounds, respectively. The trial state employed in the derivation of the upper bound can also be used in d = 2 dimensions. We refer to Proposition A.1 in Appendix A for a precise statement and its proof. A corresponding lower bound for d = 2 is still missing, however.
The analogue of Theorem 2.1 for d = 3 was proved in [13]. While the new tools developed here for the lower bound use the one-dimensional nature of the model in an essential way, they are robust enough to allow for an extension of our results to quasione-dimensional systems, like Heisenberg models defined on ladder graphs. Such an extension is rather straightforward and we shall not give the details here.

Boson representation
It is well known that the Heisenberg Hamiltonian can be rewritten in terms of bosonic creation and annihilation operators [17]. For any x ∈ [1, . . . , L] ⊂ Z, we set where a † x , a x are bosonic creation and annihilation operators, S ± x = S 1 x ± i S 2 x , and [ · ] + = max{0, · } denotes the positive part. The operators a † and a act on , and satisfy the canonical commutation relations [a, a † ] = 1. One readily checks that (3.1) defines a representation of SU (2) of spin S, and the operators S x leave the space where we denote the number of particles at site x by n x = a † x a x . It describes a system of bosons hopping on the chain [1, . . . L] with nearest neighbor attractive interactions and a hard-core condition preventing more than 2S particles to occupy the same site. Also the hopping amplitude depends on the number of particles on neighboring sites, via the square root factors in the first line in (3.2).
In the bosonic representation (3.2), the Fock space vacuum | (defined by a x | = 0 for all x) is a ground state of the Hamiltonian H L , and the excitations of the model can be described as bosonic particles in the same way as phonons in crystals. There exists a zero-energy ground state for any particle number less or equal to 2SL, in fact. While this may not be immediately apparent from the representation (3.2), it is a result of the SU (2) symmetry of the model. The total spin is maximal in the ground state, which is therefore (2SL + 1)-fold degenerate, corresponding to the different values of the 3-component of the total spin. The latter, in turn, corresponds to the total particle number (minus SL) in the bosonic language.
Before we present the proof of Theorem 2.1, we shall briefly explain the additional difficulties compared to the d = 3 case, and the reason why the proof in [13] does not extend to d = 1. Spin-wave theory predicts that at low temperatures the interaction between spin waves can be neglected to leading order. This means that (3.2) can effectively be replaced by the Hamiltonian of free bosons hopping on the lattice. At low temperature and long wavelengths 1, one can work in a continuum approximation where the last term − x n x n x+1 in (3.2) scales as −d , while the kinetic energy scales as −2 . The interaction terms can thus be expected to be negligible only for d ≥ 3, and this is indeed what was proved in [13]. This argument is in fact misleading, as the attractive interaction term turns out to be compensated by the correction terms in the kinetic energy coming from the square root factors. Making use of this cancellation will be crucial for our analysis (while it was not needed in [13] to derive the free energy asymptotics for d ≥ 3). We note that for d = 1 and d = 2 the interaction is strong enough to create bound states between magnons [15,16,20,26,27]. These occur only at nonzero total momentum, however, with binding energy much smaller than the center-of-mass kinetic energy at low energies. Hence, they do not influence the thermodynamic properties of the system at low temperature to leading order.

Upper bound
Recall the definition of C 1 in (2.3). In this section, we will prove the following.
The general structure of the proof will be similar to the corresponding upper bound given in [13]. The difference lies in the choice of the trial state, which in contrast to [13] allows for more than one particle on a single site. This is essential in order to capture the desired cancellations explained in the previous section.
Step 1. Localization in Dirichlet boxes. Our proof will rely on the Gibbs variational principle, which states that for any positive with Tr = 1. We shall confine the particles into smaller intervals, introducing Dirichlet boundary conditions. To be precise, let be the Heisenberg Hamiltonian on L : It is well known that the thermodynamic limit in (2.2) exists, hence we can assume without loss of generality that L = k( + 1) + 1 for some integers k and . By letting all spins on the boundary of the smaller intervals of side length point maximally in the negative 3-direction, we obtain the upper bound In particular, by letting k → ∞ for fixed , we have in the thermodynamic limit.
Step 2. Choice of trial state. To obtain an upper bound on f D , we can use the variational principle (4.2), with where we denote the Fock space F ≡ F for simplicity. Here, P is an operator satisfying 0 ≤ P ≤ 1, and is defined by (4.6) Note that 0 ≤ P ≤ 1, and P is zero if more than 2S particles occupy some site. The operator K is the Hamiltonian on Fock space F describing free bosons on = [1, . . . , ] with Dirichlet boundary conditions, i.e., where D denotes the Dirichlet Laplacian on and x, y means that x and y are nearest neighbors. The eigenvalues of − D are given by Step 3. Energy estimate. We shall now give a bound on the energy of the trial state.

Lemma 4.1 On the Fock space
Proof Definition (4.5) implies that It follows that With the aid of (4.10) and (4.11), one checks that The desired bound (4.9) then follows directly from P 2 1 − n x 2S 1 − n y 2S ≤ 1 and We conclude that As a next step, we will show that

Lemma 4.2 We have
Wick's rule for Gaussian states therefore implies that Moreover, Inserting this bound into (4.15) yields the desired result.
Step 4. Entropy estimate. It remains to give a lower bound on − Tr ln , the entropy of . We proceed in the same way as in [13,Lemma 4.4].

Lemma 4.3 We have
Proof We have Using the operator monotonicity of the logarithm, as well as the fact that the spectra of Pe −β K P and e −β K /2 P 2 e −β K /2 agree, we can bound Hence, (4.17) In the last term, we can bound 1−P 2 as in (4.14), and evaluate the resulting expression using Wick's rule. With φ p the eigenfunctions of the Dirichlet Laplacian, displayed below Eq. (4.8), we obtain (4.18) The expectation value of n x can be bounded independently of x as in (4.16). When summing over x, we can use the normalization x |φ p (x)| 2 = 1. To estimate the sums over p, we proceed similarly as in the proof of Lemma 4.2 to obtain In combination, this yields the desired bound.
Step 5. Final estimate. The Gibbs variational principle (4.2) together with (4.12) and Lemmas 4.3 and 4.2 implies that for a suitable constant C > 0, as long as C(β S) 1/2 ≤ (β S) 2/3 . The first term on the right side in the second line of the expression above equals By monotonicity, we can bound the sum by the corresponding integral, which is of the desired form, except for the missing part for arbitrary α > 0, some C > 0 (depending on α), and C 1 defined in (2.

Lower bound
Recall the definition (2.3) of C 1 . In this section, we shall prove the following.
Note that in contrast to the upper bound in Proposition 4.1, the lower bound above is not entirely uniform in S. Indeed, one has ln(β S 3 ) = ln(β S) + ln S 2 and hence S is not allowed to grow arbitrarily fast compared to β S. To obtain a uniform bound, one can combine our results with the method in [11] where the case S → ∞ for fixed β S was analyzed.
The remainder of this section is devoted to the proof of Proposition 5.1. For clarity, the presentation will be divided into several steps. Some of them will use results from [13].
Step 1. Localization. Recall the definition (2.1) of the Hamiltonian H L . For a lower bound, we can drop a term (S 2 − S · S +1 ) from the Hamiltonian, which leads to the subadditivity for 1 ≤ ≤ L − 1. By applying this repeatedly, one readily finds that for any ≥ 1. We shall choose large compared with the thermal wavelength, i.e., (β S) 1/2 .
Step 2. Lower bound on the Hamiltonian. Recall that the total spin operator is defined as S tot = x=1 S x . It follows from the theory of addition of angular momenta that where σ denotes the spectrum. We will use the following bound on the Hamiltonian.

Lemma 5.1 With T defined in (5.2), we have
Proof It was shown in [13, Eq. (5.6)] that for three distinct sites x, y, z, and consequently that for any x < y. After summing the above bound over all 1 ≤ x < y ≤ , we obtain We have w x=1 y=w+1 As S 2 tot = T (T + 1) we thus have The final bound (5.3) then follows from the fact that T ≤ S .
Note that Lemma 5.3 implies, in particular, a lower bound of 2S −2 on the spectral gap of H above its ground state energy. For S = 1/2, it follows from the work in [9] that the exact spectral gap equals (1 − cos(π/ )) (which is 1 2 π 2 −2 to leading order for large ).
Step 3. Preliminary lower bound on free energy. With the aid of (5.3), we shall now prove the following preliminary lower bound on the free energy. The last trace equals the number of ways n indistinguishable particles can be distributed over sites, with at most 2S particles per site. Dropping this latter constraint for an upper bound, we obtain In particular, For large β S, this expression is minimized when ≈ 0 with 0 given in (5.4). If 0 /2 ≤ ≤ 0 , we can use the lower bound on in the first term in (5.6), and the upper bound on the second, to obtain which is of the desired form. If > 0 , we can divide the interval [1, ] into smaller ones of size between 0 /2 and 0 . Using the subadditivity (5.1), we conclude (5.7) also in that case. Step Tr e −β H P E 0 ,n (5.9) where In other words, we can restrict the trace to states with S 3 tot being as small as possible (given S 2 tot ). In the particle picture discussed in Sect. 3, this amounts to particle number N = S − T = n. Because of (5.3), we have E 0 > H ≥ 2Sn/ 2 on the range of P E 0 ,n , hence the sum in (5.9) is restricted to Step 5. A Laplacian lower bound. With the aid of the Holstein-Primakoff representation (3.1), we can equivalently write the Hamiltonian H in terms of bosonic creation and annihilation operators as where n x = a † x a x ≤ 2S. Note that written in this form, the Hamiltonian H is manifestly positive, contrary to (3.2).
Let N = x n x = S + S 3 tot denote the total number of bosons. States with n particles, i.e., N = n , are naturally identified with n-boson wave functions 1 in where | denotes the vacuum (which corresponds to the state with all spins pointing maximally down). Using (5.12), we have in this representation Because of permutation-symmetry, we can also write this as For a lower bound, we can restrict the sum over x 1 , . . . , x n to values such that x k = x l for all k = . For a given j, we can further restrict to x k = x j + 1 for all k = j. In this case, the square root factors above are equal to 1. In other words, we have the lower bound where the sum is over the set X ,n := {[1, ] n : Note that we have to assume that ≥ n for the set X ,n to be nonempty. The factor 1/2 arises from the fact that particles are allowed to hop both left and right, i.e., each pair (X , Y ) appears twice in the sum. Note also that the above inequality is actually an equality for S = 1/2, since in this case no two particles can occupy the same site.
On the set {1 ≤ x 1 < x 2 < · · · < x n ≤ } ⊂ X ,n , define the map and extend it to the set X ,n = {[1, ] n : x i = x j ∀i = j} via permutations. In other words, V maps x i to x i − k i where k i denotes the number of x j with x j < x i . As a map from X ,n to [1, − n + 1] n , V is clearly surjective, but it is not injective. Points in [1, − n + 1] n with at least two coordinates equal have more than one pre-image under V . The pre-images are unique up to permutations, however, hence we can define a map V : We then have For every pair (A, B) ∈ [1, − n + 1] n with |A − B| = 1, there exists at least one pair (X , Y ) ∈ X ,n with |X − Y | = 1 in the pre-image of V . In other words, the last sum above is greater or equal to 1 if |A − B| = 1. We have thus proved the following statement. (5.13). Then,
Step 6. Bounds on the two-particle density. We will use Proposition 5.2 and the minmax principle to obtain a lower bound on the eigenvalues of H . For this purpose, we need an estimate on the norm of V .
Proof From the definition of := V , we have where |V −1 (V (X ))| denotes the number of points in the pre-image of V (X ). This number equals one if X is such that |x j − x k | ≥ 2 for all j = k. Hence, Indeed, the norm of involves a sum over all possible configurations so we need to remove the terms which correspond to x i = x j or x i = x j + 1 for some i = j. The x i = x j terms are removed through the term 1 to x i = x j + 1 are removed through −1 x=1 n x n x+1 , which is zero if and only if there are no two neighboring sites that are occupied. With = 1 and the definition of ρ (x, y) this becomes (5.14).
We shall give a lower bound on the right side of (5.14) in terms of the energy of .
The Cauchy-Schwarz inequality therefore implies that Moreover, we thus have . (5.16) We note that For given y ≤ /2, choose x y > y such that We have (where the sum is understood to be zero if x y = y + 1). The first term on the right side can be bounded as using that y ≤ /2 by assumption. For the second, we use the bound (5.16) above, which implies that |ρ (w, y) − ρ (w + 1, y)| ≤ √ 2 |h y w 1/2 (ρ (w + 1, y) + ρ (w, y)) 1/2 for w ≥ y + 1. After summing over y and w, using the Cauchy-Schwarz inequality and the fact that x,y ρ (x, y) = n(n − 1), we thus have the upper bound If y > /2, we use the symmetry of ρ and write ρ (y + 1, y) = ρ (y, y + 1) = ρ (x y , y + 1) instead, where x y is now defined by minimizing ρ (x, y + 1) for x ≤ y. Proceeding as above, we finally conclude the desired estimate.
A similar bound holds for x ρ (x, x).
Proof Since ρ (x, x) vanishes for S = 1/2, we can assume S ≥ 1 henceforth. By (5.16), It thus follows from the Cauchy-Schwarz inequality that In the last line, we can make the rough bounds 2 −1 x=1 ρ (x + 1, x) ≤ n(n − 1) and x=1 ρ (x, x) ≤ n(n − 1), and for the term in the first line we use (5.15). Using also S ≥ 1, this completes the proof of (5.17).
Step 7. Final estimate. Recall the definition (5.10) of P E 0 ,n . It follows from Proposition 5.2 that Here, we used (5.11). We shall choose the parameters such that δ 1 for large β. The min-max principle readily implies that the eigenvalues of H in the range P E 0 ,n are bounded from below by the corresponding ones of S(1 − δ)(− −n+1 n ). In particular, for any β > 0 Note that the Laplacian −n+1 n depends on n, besides the particle number, also via the size of the interval [1, − n + 1]. For a lower bound, we can increase the interval size back to , all eigenvalues are clearly decreasing under this transformation. In particular, with δ in (5.18), N 0 = E 0 2 /(2S) and E 0 = O( β −3/2 S −1/2 (ln(β S)) 1/2 ln(β S 3 )). Since ε( p) is increasing in p, we further have The error terms compared to the desired expression Looking again at (4.15), we see that the summation over x ∈ yields a factor 2 , and hence we arrive at the desired bound (A.5).
Next, we establish the two-dimensional counterpart of the entropy estimate. We have Lemma A.2 In the case d = 2, we have Tr F e −β K Tr F Pe −β K P .
Proof As in the case of the previous lemma, the only difference with regard to the one-dimensional case lies in the estimation of the p sums in (4.18). By proceeding similarly as above, we obtain ln(1 − e −β Sε( p) )d p ≤ − 1 β(β S) 1/2 ( + 1) R + ln(1 − e − p 2 )d p .