{'fraction_non_alphanumeric': 0.04476741578073835, 'fraction_numerical': 0.015655138384246176, 'mean_word_length': 4.183290707587383, 'pattern_counts': {'":': 0, '<': 42, '<?xml version=': 0, '>': 10, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Intermediate temperature dynamics of one-dimensional Heisenberg antiferromagnets

We present a general theory for the intermediate temperature (T) properties of Heisenberg antiferromagnets of spin-S ions on p-leg ladders, valid for 2Sp even or odd. Following an earlier proposal for 2Sp even (Damle and Sachdev, cond-mat/9711014), we argue that an integrable, classical, continuum model of a fixed-length, 3-vector applies over an intermediate temperature range; this range becomes very wide for moderate and large values of 2Sp. The coupling constants of the effective model are known exactly in terms of the energy gap above the ground state (for 2Sp even) or a crossover scale (for 2Sp odd). Analytic and numeric results for dynamic and transport properties are obtained, including some exact results for the spin-wave damping. Numerous quantitative predictions for neutron scattering and NMR experiments are made. A general discussion on the nature of T>0 transport in integrable systems is also presented: an exact solution of a toy model proves that diffusion can exist in integrable systems, provided proper care is taken in approaching the thermodynamic limit.


I. INTRODUCTION
One-dimensional Heisenberg antiferromagnets are strongly-interacting quantum many body systems for which a detailed quantitative confrontation between theory and experiment has been possible. A rather precise, and parameter-free understanding of their low temperature dynamic properties has emerged in a number of recent NMR experiments performed by Takigawa and collaborators [1,2,3,4]. These systems can therefore serve as useful springboards towards deciphering the behavior of interacting systems of greater complexity.
The past theoretical work on the dynamic properties of these quantum antiferromagnets has focused mainly only the universal behavior in the asymptotic low temperature (T ) regime T → 0 [5,6,7,8,4]. In the present paper we will extend the theory to a separate range of intermediate temperatures. We shall argue that under suitable conditions, to be described precisely below, this intermediate temperature range can be quite wide, and is described by a continuum dynamical model quite different from that required for T → 0. This intermediate temperature dynamics was discussed briefly for a limited class of antiferromagnets in the last section of Ref [4].
Our work will also connect with earlier investigations of the dynamics of classical lattice antiferromagnets (Refs [9,10] and references therein). In a sense, our paper provides a bridge between the modern quantum dynamics and the classical studies of the 70's. There is an overlapping window of validity for our theory and the classical investigation of Reiter and Sjölander [10], and here our results are generally consistent with theirs, although there are some details that disagree. We will review this earlier work, in the context of our results, in Section V.
A large fraction of the experimental examples of one-dimensional Heisenberg antiferromagnets consist of p parallel, coupled chains of spin S ions (for p = 1 these are ordinary spin chains, while for p > 1 these are commonly referred to as p-leg ladders). For all T < T (1) max , where T (1) max will be defined shortly, these antiferromagnets are described by a universal quantum field theory: the one-dimensional O(3) non-linear sigma model. This field theory has the quantum partition function (in units withh = k B = 1, which we use throughout) Here n(x, τ ) is a three-component unit vector representing the orientation of the antiferromagnetic order parameter at spatial position x and imaginary time τ , c is a spin-wave velocity, and H is a uniform external magnetic field-we will be interested only in the linear response to H. There are two dimensionless coupling constants in L, θ and g. The first, θ, is the coefficient of a topological term, and has the value θ = π for 2Sp an odd integer, and the spectrum of excitations above the ground state is then gapless. For 2Sp an even integer, θ = 0, and then there is gap to all excitations. The coupling g plays a role in determining the energy scale at which certain crossovers (to be discussed below) take place, but does not modify the physics otherwise. A straightforward semiclassical (large S) derivation shows that θ = 0 θ = π g 0 g g c 0 FIG. 1. Renormalization group flows for the dimensionless coupling g in (1.1). For θ = 0, g has a runaway flow to g = ∞ and the ground state is a quantum paramagnet with a gap ∆. For θ = π, there is a fixed point at g = g c , of order unity, and near it the flow is dg/dℓ ∝ (g − g c ) 2 . This fixed point is described by the k = 1, SU (2) Wess-Zumino-Witten model. The crossover between the g = 0 and g = g c fixed points takes place at an energy scale of order T 0 . The region g > g c usually corresponds to a gapped state with spin-Peierls order, and is not considered in this paper.
To reiterate, the quantum theory (1.1) applies to the lattice antiferromagnet at all T below that in (1.3).
Let us now review the well-known, T = 0, renormalization group properties of (1.1) [13]. The topological angle θ remains fixed at θ = 0, π, while the flows of the coupling g are sketched in Fig 1. For both cases θ = 0, π there is a fixed point g = 0 which is unstable at low energies. Indeed the beta function describing the flow away from g = 0 is independent of θ to all orders in g. However, non-perturbative, topological effects do distinguish the two values of θ. For θ = 0, the flow is believed to continue all the way to g = ∞, corresponding to a quantum paramagnetic ground state with an energy gap ∆. In contrast, for θ = π, the flow is into a strong coupling, infrared stable fixed point at g = g c . There is a scale-invariant and gapless theory which describes this fixed point-the k = 1, SU(2), Wess-Zumino-Witten model. For g > g c , there is a runaway flow to g = ∞, usually associated with the appearance of spin-Peierls order; this last regime will not be discussed in this paper.
Our primary interest here shall be in the region in the vicinity of the unstable g = 0 fixed point. For both θ = 0, π, there is a characteristic energy scale, usually denoted Λ M S , which determines the location of the crossover from the vicinity of the g = 0 fixed point to the strong-coupling behavior. For energies or temperatures smaller than Λ M S , the strong coupling behavior should apply, and as we have just discussed, this is quite different for θ = 0 and θ = π. However, the physics at energies or temperatures larger than Λ M S is controlled by the flow in the vicinity of the g = 0 fixed point, and this is common to both θ = 0 and θ = π. The energy scale Λ M S can be estimated from the structure of the perturbative β-function. For small g we have neglected here a prefactor of a power of g coming from higher-loop corrections. So from (1.2), for either S or p moderately large, the scale Λ M S becomes exponentially small.
We are now ready to discuss the static thermodynamic properties of the quantum field theory (1.1) as a function of T . We will characterize the system by T dependence of two important observables: χ u (T ) and ξ(T ). The first is the uniform susceptibility, χ u , which is the linear response to the field H = (0, 0, H): χ u = (T /L)(d 2 ln Z Q /dH 2 )| H=0 , where L is the (infinite) length of the spatial direction; this is the susceptibility per rung of the ladder. The second is the correlation length, ξ(T ), which determines the exponential decay of the equal-time two-point n field correlator as a function of x.
We will consider temperatures above and below Λ M S in turn. (A) T < Λ M S For T < Λ M S , as just noted, we must distinguish θ = 0 and θ = π. (A.I)θ = 0 For θ = 0, there is an energy gap ∆, and the susceptibility is simply that of a dilute, thermally activated, classical gas of triplet magnons above the gap; these contribute an exponentially small susceptibility [14,15] Experimentally, we can view (1.5) as the definition of the gap ∆ and the velocity c, which are to be determined by fitting measurements to (1.5). The correlation length, ξ(T ) takes a finite, T -independent value in this quantum paramagnet, up to correlations exponentially small in ∆/T ; for the case where ∆ is significantly smaller than J, we have [12,16] For the gapless case, θ = π, there are excitations with non-zero spin at arbitrarily low energies, and so χ u (T ) remains non-zero as T → 0 [18,19]: Again, this experimentally defines c and a new temperature scale T 0 which determines the onset of a logarithmic correction to the T = 0 susceptibility due to the slow flow into the fixed point at g = g c . The correlation length of this critical paramagnet now diverges as T → 0 [5,6,7,8,20]: The dynamical properties of quantum antiferromagnets in the low temperature regime T < Λ M S have been discussed at length in Ref [4] for the gapped case (θ = 0) and in Refs [5,6,7,8] for the gapless case (θ = π).
Note that both ∆ and T 0 are energy scales characterizing the flow into the strong-coupling region. These should therefore be universally related to Λ M S which is the scale of flows out of the weak coupling region. We will discuss the universal relation shortly, once we have defined Λ M S more precisely. (B)T > Λ M S Let us now consider the regime Λ M S < T < T (1) max ; the upper-bound is necessary to ensure the continuum quantum theory still applies. The existence of this intermediate temperature regime requires that Λ M S < T (1) max , a condition that is not well satisfied for small S and p, and so this regime almost certainly does not exist for p = 1 and S = 1/2, 1, but there is evidence that it is present for p = 1, S = 2 [17]. Here we are controlled by physics in the vicinity of the g = 0 fixed point, and it should be possible to treat quantum fluctuations in a renormalized perturbation theory in g. As discussed in Ref [4], a non-perturbative treatment of the thermal fluctuations is still necessary, but this can be carried out exactly because of the low spatial dimensionality. The result of such a calculation is [4] where γ is Euler's constant. This result can also be viewed as the precise experimental definition of Λ M S . For completeness, we also quote the result in this regime for the correlation length, ξ(T ): (1.10) We now give the promised universal relationship between the energy scales characterizing the weak (Λ M S ) and strong (∆, T 0 ) coupling regimes. This requires non-perturbative knowledge of the renormalization group flows, and can only be obtained from an analysis of the full thermodynamic Bethe-ansatz solution of the quantum field theory (1.1). For θ = 0, such an analysis was carried out in Refs [21], and the result is now well known: For the gapless case, θ = π, we will present a derivation of the required relationship in Appendix A, building upon some recent results [23,22]; our result is: (1.12) We are now finally in a position to state precisely the main objective of this paper. We will describe the dynamical properties of one-dimensional antiferromagnets in the intermediate temperature regime Λ M S < T < T (1) max . We quickly note as an aside that most of our results actually hold over a wider regime of temperatures Λ M S < T < T (2) max , where we will define the discuss the origin of T (2) max > T (1) max below; for now we ignore this point. A common treatment is possible for the θ = 0 and θ = π in this regime, with the two cases differing only in the input values of the static parameters χ u (T ) and ξ(T ) as defined by (1.9,1.10,1.11,1.12). Moreover, as Λ M S becomes exponentially small for moderate values of S or p, this regime can be quite wide and should be readily observable experimentally. Indeed, there is good evidence from recent measurements of static properties in quantum Monte Carlo simulations [17] that this universal intermediate temperature regime exists even for S = 2 spin chains.
The formulation of the dynamics properties for Λ M S < T < T (1) max was already discussed in Ref. [4]. The key point [4] is to notice that the energy of a typical spin-wave excitation, which is of order cξ −1 (T ), is parametrically smaller than T when ξ obeys (1.10). So the thermal occupation number of these spin-wave modes is large: The second expression in (1.13) is the classical equipartition value, which indicates that the spin-wave excitations may be treated classically. The classical partition function controlling these fluctuations can be deduced by demanding that its correlations match with (1.9,1.10), while the dynamic equation of motions following by replacing the quantum commutators with Poisson brackets. In this manner, the problem reduces to the effective classical phasespace partition function [4] (1.14) Here n(x) is a classical variable representing the orientation of the antiferromagnetic order and L(x) is its classical, canonically conjugate angular momentum. Because n(x) is of unit length, its motion is always in a direction orthogonal to its instantaneous direction, and there is no radial kinetic energy: the square of angular momentum represents the entire kinetic energy, and χ u⊥ is the moment of inertia of the fluctuating n. The value of χ u,⊥ can be determined by realizing [4] that L(x) is simply the classical limit of the quantum operator corresponding to the magnetization density then demanding that Z C reproduce the correct uniform susceptibility to a uniform external field under which H C → H C − dx H · L, we obtain the factor 2/3 comes from the constraint L · n = 0, so there are only two independent components of L at each spatial point. Notice that (1.14) involves a functional integral over the commuting fields n and its conjugate momentum L as functions only of the spatial co-ordinate x, but is independent of the real time, t. It therefore yields only equal-time correlation functions, as is the usual situation in classical statistical mechanics. To obtain unequal time correlators, we have to separately specify the equations of motion, and these are obtained by replacing quantum commutators with Poisson brackets. For the fields n(x), L(x) these are where α, βγ = 1, 2, 3. The equations of motion (in real time) now follow from the Hamiltonian H C , and they are We are now interested in unequal time correlation functions of (1.18), averaged over the classical ensemble of initial conditions specified by Z C .
To complete the quantum-to-classical mapping, we recall [4] the relationship between the correlations of the underlying antiferromagnet and the quantum field theory Z Q , and those of the classical non-linear wave problem defined by (1.14) and (1.18). Correlations of the antiferromagnet in the vicinity of the antiferromagnetic wavevector are given by the correlations of n(x, t) under Z Q , and are related to those in the classical problem by [4] n(x, t) · n(0, 0) where the subscript Q represents averages under the quantum partition function Z Q , the subscript C represents averages under the classical dynamical problem defined by (1.14,1.18), and A is an overall T -independent normalization related to the amplitude of the correlations at T = 0. Next, correlations of the antiferromagnet in the vicinity of zero wavevector are given by the correlations of L(x, t) under Z Q , and are essentially equal to those in the classical problem by [4] L(x, t) · L(0, 0) Q = L(x, t) · L(0, 0) C . (1.20) The absence of an overall rescaling factor here is related to the conservation of the total magnetization density.
The main objective of this paper is to evaluate n(x, t)n(0, 0) C and L(x, t)L(0, 0) C . An important property of these correlators is that they satisfy simple scaling laws which allow us to completely scale away all dependencies on ξ(T ) and χ u⊥ (T ), and to express everything in terms of parameter-free, universal functions. These scaling laws follow from the fact that (1.14,1.18) define a continuum classical problem which is free of all ultraviolet divergences: this will become evident from our analytic computations in Section II and the numerical results of Section III. Consequently, simple engineering dimensional analysis involving rescaling of x, t and L can be used to absorb dependencies on the dimensionful parameters. In this manner, it is not difficult to show that [4] n(x, t) · n(0, 0) C = Φ n (x, t) and Φ n (x, t) and Φ L (x, t) are universal scaling functions. It is the primary task of this paper to determine these scaling functions. In principle, these scaling functions are determined by solving (1.14,1.18) after setting all parameters equal to unity, T = ξ = χ u⊥ = 1, while replacing x,t by x,t. Notice that the resulting equation is then parameter-free, and so there is no explicit small parameter in which any kind of expansion can be carried out; nothing short of an exact solution will do. We think this reasoning invalidates some of the conjectures on exactness of results made in Ref [10], as we will discuss further in Section V. However, when the arguments of the scaling functions are themselves small, i.e. |x|, |t| ≪ 1, then a systematic perturbation expansion is possible, and this will be presented in Section II. We close this introductory discussion by returning to the issue of the maximum temperature up to which these results can be applied to lattice antiferromagnets. The appearance here of a classical spin model suggests that one should think about classical spin models obtained by starting directly from the lattice quantum spin model, without the use of the quantum field theory Z Q as an intermediate step [11]. Such a classical description will only work for large S, and we can ask the question of when the resulting classical model can be described by a continuum classical theory, which will clearly be the one defined by (1.14,1.18). The correlation length of a classical spin antiferromagnet is of order JS 2 pa/T , and a continuum description will work provided this is larger than the lattice spacing a and the width of a ladder system pa. This gives us the estimate [11] T (2) max ∼ JS 2 .
(1. 23) In the regime T (1) max < T < T (2) max , we can use a purely classical description of the spin model: its correlation length and uniform susceptibility will not be universal, but has to computed for the specific model under consideration. For the model with only nearest neighbor exchange, a standard computation on the classical antiferromagnet gives for the uniform susceptibility and for the correlation length, It is satisfying to note that there is a precise agreement between the results (1.9,1.10) and (1.24,1.25) at the common boundary of their respective regions of validity, T ∼ T (1) max : for the nearest neighbor model under consideration we use the estimates for g and c in (1.2), and then using ln(T /Λ M S ) ≈ 2π/g, we find the required agreement.
It appears useful to review the final status of the regimes of validity of the model studied here. The universal, continuum, classical model (1.14,1.18) describes all one-dimensional Heisenberg antiferromagnets in the temperature regime Λ M S < T < T (2) max ; this regime is wide and well defined for moderately large values of S or p, and there is evidence that it exists already for p = 1, S = 2 [17]. The complete definition of this classical dynamical model requires the input of the temperature-dependent static parameters ξ(T ) and χ u (T ). In the regime Λ M S < T < T (1) max these parameters are universally specified by (1.9,1.10), with Λ M S given by (1.11) for gapped spin chains (θ = 0) and by (1.12) for gapless spin chains (θ = π). In the higher temperature regime T (1) max < T < T (2) max , these parameters are given by (1.24) and (1.25) for the model with only nearest neighbor exchange, and by related non-universal expressions for other antiferromagnets.
The following sections contain technical details towards the determination of the scaling functions Φ n and Φ L in (1.21), along with some theoretical analysis on the relationship between integrability and diffusion. We will begin in Section II by describing the analytical short time expansion of the correlators Φ n and Φ L . The long time limit will then be studied numerically in Section III. The subsequent Section IV discusses issues which are somewhat peripheral to the main focus of this paper: the continuum equations of motion (1.18) are known to be integrable [24], and this raises numerous fundamental questions on the nature of spin transport in integrable systems. These will be addressed in Section IV by the study of a simple, integrable, toy model whose spin correlators can be determined in close form. Further, we will see that these correlators have a striking similarity to those of (1.18). Readers interested primarily in spin chains can omit Section IV and skip ahead to Section V where we will describe the implications of our results for experiments.

II. SHORT TIME EXPANSION
This section will determine the small t expansion of the scaling functions Φ n and Φ L in (1.21). Normally, there is a completely straightforward way of determining the short time expansion of interacting systems [25]-it can be related, order-by-order, to equal-time correlators involving higher moments of the fields. However, this standard procedure does not work for the model (1.14,1.18) of interest here. This is because we are dealing with a continuum model with an infinite number of degrees of freedom, n(x), L(x), present at arbitrary short distance scales. If we naively generate the moment expansion, we find that the terms quickly acquire rather severe ultraviolet divergences.
A separate theoretical tool is necessary to generate the short time expansion, and this shall be described here. We shall use an analog of the field-theoretic method known as chiral perturbation theory. As we shall see below, the expansion is actually in powers of |t|-this implies a non-analyticity in the t-dependence at t = 0, which is in fact the reason for the ultraviolet divergences in the moment expansion. The latter method only gives an analytic expansion in t, by construction.
In this section, and in Appendix B, we will use units in which T = χ u⊥ (T ) = 1. However, we will retain explicit dependence on ξ ≡ ξ(T ). It turns out to be quite useful to keep track of powers of ξ. Indeed, our computations will be designed to generate an expansion of the correlators in powers of 1/ξ, and this is a posteriori seen to be a short time expansion. We will return to physical units in stating our final results.
We will therefore consider the problem of unequal time correlation functions of the nonlinear partial differential equations when averaged over the ensemble of initial conditions defined by the partition function The last term in the action represents a field of strength ξm 2 turned on in the z direction and serves as a regularization parameter for our perturbation expansion. At the end of our calculations we shall let m → 0. From this partition function, one can immediately find the equal time correlation functions to be the subscript C is implied on all averages in this section, unless stated otherwise. We first construct our perturbation expansion in powers of 1/ξ for the equal time problem specified by (2.2), and check that we do arrive at the correct correlation functions as specified in (2.3). The extension to the unequal time problem will then be straightforward. First, the constraints on the fields n and L are solved by introducing two complex scalar fields φ and ψ.
We introduce this decomposition to the functional integral (2.2) and expand the square roots in power series of 1/ξ. In this manner, we arrive at an interacting field theory with an infinite number of interactions, with 1/ξ as the small coupling. However, to any particular order in perturbation theory in 1/ξ, we only need to keep a finite set of interactions terms. We evaluate the correlation functions in real space using the ordinary machinery of diagrammatic perturbation theory. It is well known that such a perturbation expansion is plagued by infrared divergences. The decomposition (2.4) takes it for granted that the n field is ordered in the z direction and φ represents small spin wave fluctuations around the ordered state. The infrared divergences are a signature of the fact that this assumption is wrong in one dimension. By introducing the external field ξm 2 , we introduce long range order into the system and thus regularize the divergences. The divergences show up as poles in 1/m. But if we calculate the O(3) invariant correlation functions as in (2.3), we find that all poles cancel. Keeping this in mind, it is not difficult to show that the results of perturbation theory agree with (2.3) at every order. The details are to be found in Appendix B.
Once we have the correlation functions in the equal time ensemble, we proceed to evaluate the unequal time correlation functions. To do that we again insert (2.4) into the equations of motion (2.1) and expand in powers of 1/ξ to get Higher order terms add further nonlinear interactions. We solve the initial value problem for each of the fields by an iterative strategy. First the free wave equation is solved and the solution is plugged back into the lowest order nonlinear term to solve the problem to the first order. To evaluate the correlation functions, we just multiply the fields and carry out the average over initial conditions. The initial condition averages are, of course, known from the calculations described above. The technical details are relegated to Appendix B and here we only quote the results to one loop order: where x and t are defined in (1.22), and the delta function in In these results if we set t = 0, we immediately recover the equal time results (2.3) to the corresponding order in 1/ξ. The structure of the correlation function (2.6) reflects the causal propagation of the conserved angular momentum L. The first term simply represents the free propagation of the angular momentum density which is completely concentrated on the "light cone". However, interactions to order 1/ξ do modify the free propagation and transfers the angular momentum density from the surface of the light cone to its interior. This is represented by the second term which is non-vanishing only inside the light cone. It can be checked that the spatial integral of (2.6) remains independent of t, as must be the case due to conservation of total angular momentum.
The result (2.7) has significant implications for the dynamic structure factor, S(k, ω), of the antiferromagnetic order parameter: (2.8) First, we note that the result (1.19) and the scaling form (1.21), imply that the dynamic structure factor satisfies where S(k) is the equal time structure factor, which is known exactly (apart from the overall normalization A) and Φ S is a universal scaling function of which describes the relaxation of the equal time correlations. The prefactor in (2.9) has been chosen so that the frequency integral over Φ S is normalized to unity for every k We will now show that, modulo some very mild assumptions, the result (2.7) exactly fixes the form of Φ S (k, ω) for |k| ≫ 1. In this short distance regime, we are at distances shorter than the correlation length, and the system should look almost ordered. So we may expect that the spectrum consists of weakly damped spin waves, and this motivates the following ansatz for Φ S in the regime |k| ≫ 1 and |ω| ∼ |k|: where γ(k) is the unknown spin-wave damping parameter; we will shortly determine the large k limit of γ(k). We now have to take the Fourier transform of (2.9), (2.10) and (2.13), compare the result with (2.7). In making this comparison, we should keep in mind that the large k, ω behavior of S(k, ω) can only determine the non-analytic terms in the small x, t expansion. First, the integral over frequencies can be performed exactly for the form (2.13), and (2.9) implies We reiterate that (2.13) is valid only for |k| ≫ 1, and so only the 1/k 2 contribution from the (k 2 + 1) denominator in (2.14) can be taken seriously. Also, we already have a non-analytic |t| dependence in the exponential, and so we can expand this in powers of |t|; in this manner we reduce (2.14) to where k 1,2 are some large positive constants, and we have assumed that γ is an even function of k. The values of the integrals over k surely depend upon p and p ′ , but the key observation is that the non-analytic terms in x and t do not. This follows from the result where all omitted terms can be written as a series in non-negative, even integer powers of x (this allows an additive constant, independent of x). We now assume γ(k → ∞) ∼ k α , and then demand consistency between (2.15,2.16) and the non-analytic terms in (2.7). It is easy to see that we must have α = 0, and so γ(∞) is a constant. Further, the unknown additive constant associated with the second integral in (2.15) must be such that there is no single |t| term in the correlator. Applying (2.16) to (2.15) with this understanding, we deduce that where all omitted terms are either analytic in x and t, or involve subleading non-analyticities. We should now compare (2.17) with (2.7) by matching only the non-analytic terms in the vicinity of the light cone x = ±t. In this latter region we can approximate |x + t||x − t| in (2.7) by 2|t||x ± t|-then the non-analytic terms in (2.17) and (2.7) match perfectly, and we obtain one of our important exact results Returning to physical units via (2.11), we conclude that the frequency Γ(T ) ≡ (ω/ω)γ(∞), given by describes damping of spin-waves for |k|ξ(T ) ≫ 1. We will discuss the experimental implications of this result in Section V.
Results for the spin-wave damping have been obtained earlier by Reiter and Sjölander [10], in their studies of classical lattice antiferromagnets. If we insert the classical values (1.24,1.25) into (2.19), we obtain which agrees with their result for p = 1. Keep in mind, though, that our result (2.19) has a much wider regime of applicability, beyond temperatures in which a purely classical thermodynamics holds. Further discussion on the relationship between our and earlier results appears in Section V.

III. NUMERICAL RESULTS
The previous section allowed us to determine the two-point correlators for small |t|. Here we will present numerical simulations which examine the large |t| limit. These were performed on a discrete lattice realization of (1.14,1.18), with lattice spacings of ξ/16 or larger, and had only nearest neighbor couplings between the n vectors. Initial states were generated by thermalizing the system by the Wolff algorithm [26]. The time evolution was carried out by a fourth order predictor-corrector method, and its accuracy was tested by keeping track of the conserved total energy and the lengths of the n vectors.
Our simulations are similar to many earlier studies of classical spin chains (see Refs [27,28] and references therein). However, there is an important difference in that we are dealing with rotor variables n, L, rather than classical spins S which obey Poisson bracket relations like those for L in (1.17).
We will consider correlators of n and L in the following two subsections.

A. Correlations of n
We first consider dynamic correlators of n. The aim of our simulations is to obtain results for the dynamic structure factor, S(k, ω), in regimes beyond the case |k| ≫ 1, |ω| ∼ |k| which was studied by the short time expansion. Our results were obtained for two cases-at equal positions (local) and at zero wavevector.
The local correlator is measured in NMR experiments, and we computed the local dynamic structure factor, S l (ω) defined by where Φ l is a fully universal function (with no arbitrariness in its overall amplitude or the scale of its argument) with a unit integral over frequency Our results for the local time-dependent correlations are shown in Fig 2 and its Fourier transform to frequency in Fig 3. Two different lattice spacings were used, and the good overlap of the data confirms that we are examining the continuum limit. The correlations decay rapidly in time, but also show a brief, but clear oscillation; this oscillation results in a finite frequency peak in S l (ω). We will discuss the physical origin of this oscillation after we have considered the zero momentum correlator.
Turning to the zero momentum correlator, we express our results in the scaling form (2.9) and obtain values for the scaling function Φ S (0, ω). We emphasize that the ansatz (2.13) does not hold for k = 0. Our results for the k = 0 correlators of n are shown in the time domain in Fig 4 and after the Fourier transform to frequency in Fig 5. As with the x = 0 correlations above, the k = 0 correlator shows a rapid decay, along with a brief oscillation; the latter leads to a finite frequency shoulder in Φ S (0, ω).
How do we understand this finite oscillation frequency observed in both the x = 0 and k = 0 correlators of n ? One way is to compare with the exactly known results [29,30] of the model with an N-component vector n in the limit of large N. At N = ∞, S(0, ω) consists of a delta function at a finite frequency ω ∼ T / ln(T /Λ M S ). So we can view the finite frequency as a remnant of the N = ∞ response at N = 3. However, there is a related, more physical, way to interpret it. The underlying degree of freedom have a fixed amplitude, with |n| = 1. However, correlations of n decay exponentially on a length scale ξ(T )-so if we imagine coarse-graining out to ξ(T ), it is reasonable to expect significant amplitude fluctuations in the coarse-grained field, which we call φ α . On a length scale of order ξ(T ), we expect the effective potential controlling fluctuations of φ α to have minimum at a non-zero value of |φ α |, but to also allow fluctuations in |φ α | about this minimum. The finite frequency in Figs 3, 5 is due to the harmonic oscillations of φ α about this potential minimum, while the dominant peak at ω = 0 is due to angular fluctuations along the zero energy contour in the effective potential. This is interpretation is also consistent with the large N limit, in which we freely integrate over all components of n, and so angular and amplitude fluctuations are not distinguished. . We used a lattice spacing ξ/8. By (2.10), the t = 0 value of this should be 2, and the difference is due to the finite lattice spacing.

B. Correlations of L
We obtained numerical results only for the x = 0 correlator of L. The short time behavior of this is given in (2.6). At long times, we expect the conservation of total L to be crucial in determining its asymptotic form. In particular, one natural assumption is that the long time correlators of L are diffusive; in this case we expect at large t. Consistency of this with the scaling form (1.21), implies that the diffusion constant, D must obey where B is a dimensionless universal number. Our numerical analysis of the autocorrelation was carried out on a system of 800 sites. The predictor-corrector method turns out to exactly conserve angular momentum and we maintained energy conservation to 4 significant digits over the duration of the simulation. We averaged over 9,600 initial conditions.
First we tested our results against the known exact short time expansion. This is shown in Numerical results for (ξ(T )/T χ u⊥ (T )) L(0, t) · L(0, 0) as a function of t in a log-log plot. The smoother line is for lattice spacing ξ/8, while the noisier line is for lattice spacing ξ/16. The agreement of the two results is evidence that this decay is a property of the continuum limit. A straight line fit (not shown) to the ξ/8 data is almost perfect, and its slope indicates that the correlator decays as t −0.61 . An equally good fit to the data was obtained by the function a 1 / √ t + a 2 /t, with the second subleading term contributing only about a 10% correction at the largest t.
L(x, t) · L(0, 0) = 2 where ω k = (2(1 − cos(kǫ))/ǫ 2 ) 1/2 . It can be verified that (3.6) reduces to (2.6) in the limit ǫ → 0. As is clear from Fig 6, the agreement between the analytical and numerical computations is quite satisfactory. Finally, we turn to the numerical results at large t. These are shown in Fig 7. A best fit to the data with a power law t −α gave an optimum value of α = 0.61. However an equally good fit to the data was obtained by the function a 1 / √ t + a 2 /t, with the second subleading term contributing only about a 10% correction at the largest t (we found a 2 /a 1 = 0.65). This second fit is consistent with diffusion-assuming this is the correct form, we obtain the estimate for the numerical prefactor in (3.5): In their studies of the classical lattice model, Reiter and Sjölander [10], also computed the spin diffusivity. Diffusion is a property of the |k| ≪ 1 regime, and we do not expect their perturbative techniques to be exactly valid. Combining the classical values of χ u (T ), ξ(T ), in (1.24,1.25) with (3.5), their result translates into the value B = 1. This value is clearly inconsistent with our numerical result above.

IV. INTEGRABILITY AND DIFFUSION
This section will examine an integrable toy model of spin transport; readers interested mainly in the experimental implications of our results so far can move ahead to Section V. Others, not interested in the details of the toy model, may want to jump to Section IV C where we will discuss general implications of the toy model solution on spin transport in integrable systems.
The toy model we shall introduce is a variant of an effective model, considered in Ref [4], for the dynamics in the regime T < ∆ for gapped chains. Here, our strategy will be to introduce the model as worthy of study in its own right, as it is simple enough to allow determination of the spin density correlator in closed form at all times. In the long time limit, the correlator has a diffusive form, and so this example proves that there is no general incompatibility between integrability and diffusion. Further motivations in examining this model are: (i ) We shall show that the short time behavior of the toy model is very similar to our result (2.6) for the continuum wave model (1.14,1.18). This is suggestive, and indicates that the long time diffusive behavior in (3.4) is not an unreasonable postulate. (ii ) The correlators of the integrable toy model can also be studied for the case of finite system of size L with periodic boundary conditions. This allows us to carefully examine the interplay of the limits t → ∞ and L → ∞. The diffusive form only appears if the L → ∞ is taken first. In the opposite order of limits we find recurrent behavior with a great deal of structure dependent upon the microscopic details of the model. We think this issue of the orders of limit of t → ∞ and L → ∞ is of considerable relevance to recent studies of T > 0 transport in integrable systems [31,32,33,34,35], and this will be discussed further in Section IV C.
We begin by describing the toy model. Place N point particles of equal mass at positions x i (i = 0, 1, . . . N −1) which are chosen independently from a uniform distribution on a circle of length L. Now independently give each particle a velocity v i , drawn from some distribution g(v), and a 'spin' m i , drawn from some distribution h(m). As the system evolves, the particles will move in straight trajectories, transporting their spin along with them. This will happen, until two particles collide, and we now have to describe the nature of such collisions. We will restrict the collisions to satisfy the important constraints of conservation of total energy, momentum, and spin in each collision. The first two are already sufficient to determine the fate of the velocities (see Fig 8): if we consistently label the particles from left to right (i.e. as we move around the circle anti-clockwise, we always encounter the particles in the order x 0 , x 1 , x 2 , . . . x N −1 ), then the particles will simply exchange velocities in each collision. In other words, in a collision between particle i and particle i + 1, the velocity of the particle i after the collision is that of particle i + 1 before the collision, and vice versa. How about the fate of the spins m i and m i+1 ? In principle, we can choose numerous possibilities interpolating between zero to total reflection, consistent with conservation of total spin i m i . The exactly solvable models are the two extremes: zero or total reflection. The case of zero reflection is rather trivial and leads only to simple ballistic transport of spin along straight lines. We will therefore consider only the case of total reflection here: in the convention we are following of labeling the particles here, this corresponds to the statement that each m i is a constant of the motion (see Fig 8). So to summarize: in each collision the particles exchange velocities but their spins 'bounce off' each other. We note that, although no exact solution is known, we expect the long-time correlations of a model with only partial reflection to be quite similar to that of the total reflection case, but with renormalized transport coefficients.
We shall be interested here in computing the correlators of the 'spin' density, L(x, t), defined by where x i (t) are the positions of the particles: these consist of piecewise straight lines which reflect at each collision. We shall compute the two-point correlator of L(x, t), averaged over the ensemble of initial conditions defined above. Further, we will choose our distributions g(v) and h(m) to be even i.e. on the average, the total net momentum and spin are zero. Then, because the initial momenta and spin are uncorrelated, we have where P (x, t) is the probability that a particle at x = 0 at time t = 0 is at the position x at time t, and ρ = N/L is the density of particles. Most of our results will be on a particular, simple, velocity distribution, which is designed to mimic the properties of the continuum model (1.14,1.18): (4.4) So each particle is allowed to have only one of two velocities, ±c. This is similar to the fact that linear spin-waves in the continuum wave model also have velocities ±c. The remainder of this section will describe the computation of the function P (x, t) using the method of Jepsen [36]. Although later more elegant solutions were put forward by Lebowitz and Percus [37] for solving the model in the thermodynamic limit, we shall use the relatively cumbersome machinery of Jepsen because it allows us to consider finite systems.
At time t = 0 the N particles are at random positions on the ring. We shall put the origin of the coordinate system at the location of particle 0. The rest of the particles from 0, 1, 2, · · · N − 1 are numbered such that the (i + 1)-th particle is immediately to the right of the i-th one. As a particle moves with uniform velocity, in a "space-time" diagram we can represent its motion as a straight line which we'll call a trajectory. When two trajectories cross, there is a collision. The particles bounce off each other, and in effect the particles exchange trajectories. So at the beginning the zeroth particle starts on the zeroth trajectory and as this trajectory crosses others, the zeroth particle moves onto a different trajectory. Now define A jk (t) to be one if the particle j is on the trajectory k at time t and zero otherwise. If we form an ensemble of systems, the the average A jk (t) is the probability that the particle j is on the trajectory k at time t. A knowledge of A jk (t) for all values of j, k, t constitute a full solution to the dynamics of the system. The solution is defined by [36]  Also, for the distribution (4.4), noting the fact that |v k − v h | = 0, 2c, we arrive at Let us define a time T ≡ L/c which is the time required by a free particle to go once around the system. Every trajectory returns exactly to its starting point after this interval of time.
It is now obvious that in a period of time NL/c = NT , each particle will return to its initial positions and velocities. Thus the Poincare recurrence time of the system, with the velocities chosen under (4.4), is of order N, rather than being of the order of e N or larger.

A. Diffusion in the Thermodynamic Limit
Now we go to the thermodynamic limit and address the question of diffusion. The limit is defined such that N and L approach infinity while the density N/L = ρ being held finite.
The zeroth particle starts out at the origin and we ask what is the probability, P (y, t), that it is at position y at time t (both y and t are held finite as the limit L → ∞ is taken). This can be written as where the angular brackets an average over all possible initial ensembles of velocity and positions of particles, while keeping the zeroth particle at the origin. The average can be evaluated exactly, and there is a simple, closed form result: . (4.9) I 0 and I 1 are the modified Bessel functions of order zero and one respectively. The resemblance to the correlation function of the non-linear wave model given in (2.6) is clear. The first term in (4.9) is a delta function along the light cone, but its contribution decreases exponentially with time. The second term lies within the light cone, and becomes increasingly important for large time. Also, if we take the short time limit of (4.9), we get which is precisely of the form (2.6). However, unlike (2.6), we can now also study the long time limit analytically. We take this limit within the light cone with y ∼ √ t, and then the asymptotic expansions of the modified Bessel functions yields which is the diffusive form assumed for the classical wave model in (3.4). As shown by Jepsen, this calculation can also be done for a general velocity distribution g(v), and provided the distribution is symmetric in v, we obtain (4.12) but with

B. Effect of Periodic Boundary Conditions in a Finite Geometry
Here we go back to the finite system and think more about it. The fact that the Poincare recurrence time is only linear in N, is because of the fact that the phase space becomes very restricted once we allow only two possible velocities. An interesting effect is that the recurrence time can be very different depending upon whether there are even or odd number of particles in the system. The recurrence time is the lowest common multiple of T (the time required for trajectories to return to their initial position) and NL/2c = NT /2 (the time required for particles to come back to their initial trajectory). With an odd number of particles, N = 2p − 1, the recurrence time is NT = (2p − 1)T as mentioned above. But if we add one single extra particle to the system, the recurrence time almost becomes half because NL/2c = NT /2 = pT is an exact multiple of T .
In fact, the recurrence time could be even smaller. If we take an even number N of particles and choose their velocities as ±c randomly, the most likely scenario is one where half of them have velocity +c, and the rest have velocities −c. In that case Amongst all the factors that we have in the product on RHS, half of them will not contribute any phase to the product because for them v k − v h = 0. The other half will contribute a phase of exactly exp(±2iu) each. So the total phase contribution will be exp(±2iuN/2) ≡ 1! Hence the most probable recurrence time for a random ensemble of even number of particles is T . It turns out that with reflecting hard wall boundary conditions, the recurrence time is 2T regardless of the initial conditions. Now let us see how the effect of the recurrence time might show up in the probability distribution of a diffusing particle. We shall essentially try to find the autocorrelation function for a single specific particle moving around in a ring for all times t < T = L/c and t > T = L/c. We shall again choose this particle to be the zeroth particle and at time t = 0, its position and velocity are x 0 ≡ 0, v 0 respectively. The probability distribution is defined as before in (4.8). Since A jk (t) is periodic with a period of NT /2, all distribution functions will also be periodic with the same period. Moreover When we carry out the average over the initial conditions, the particles k and h will have the velocities +v k , +v h and −v k , −v h with equal probability. So the probability distribution function will satisfy P (y, NT /2 − t) = P (y, t). Thus we need to evaluate this function only for 0 ≤ t ≤ NT /4. The details of the evaluation of P (y, t) are again relegated to the appendix. Let us represent the distribution function as P (y, t) = [δ(y + ct) + δ(y − ct)]P (1) (t) + [θ(y + ct) − θ(y − ct)]P (2) (y, t). (4.17) We write down explicit series solutions for P (1) (t) and P (2) (0, t). These sums could not be evaluated in a closed analytic form. So we carried out the sums numerically for specific  values of N and L. P (1) (t) and P (2) (0, t) are plotted in Figs 9 and 10 respectively. For times t ≪ T and N ≫ 1 with N/L fixed, these results reduce to the results (4.9) derived in the thermodynamic limit. For t > T , we see lots of complicated structures. Not all of them are well understood. But some of the prominent ones are easy to understand. For example at times t = nT , n = 1, 2, 3, . . . we see peaks in P (1) , each of whose height turns out to be 1/ √ 2πN. The origin of this is very easy to understand. Out of the whole ensemble of initial conditions, a fraction of them will have exactly half of the particles with velocity +c and the others with velocity −c (with the assumption that N is even). Since these set of initial conditions have a recurrence time of only T , so they come back to the original distribution after a time T and hence contribute to this peak. Using the binomial distribution, it is easy to prove that the fraction of ensembles which have exactly half the particles with one velocity and the other half with the opposite velocity is exactly 2/πN. The extra factor of 1/2 comes from the fact that we have both left and right going initial conditions for our test particle. If we change the number of particles by one, N becomes odd. So no set of initial conditions will have a recurrence time T , and these peaks will disappear as shown in Fig 11.

C. Discussion
Our study of an integrable toy model in this section has highlighted the extreme importance of taking the limits of long time t → ∞, and large system size L → ∞ with proper care.
If we send L → ∞ first, then we explicitly demonstrated the existence of spin diffusion in the subsequent long time limit. This is, of course, the correct thermodynamic limit, and the presence of spin diffusion also makes physical sense: once the limit L → ∞ has been taken, an infinite number of parameters are needed to specify the initial thermal state of the system-diffusion then arises when some local degree of freedom starts to sample an increasing number of these infinite number of random initial conditions with the passage of time.
On the other hand, very different results were obtained for the limit t → ∞ at any fixed L. Here the integrable nature of the system was immediately evident, and we observed a bizarre set of recurrences dependent sensitively on details of the microscopic Hamiltonian. The Poincare recurrence time of our toy model was quite short, and this was clearly due to its integrability. No sign of spin diffusion was seen.
A number of recent studies have examined the issue of T > 0 spin transport in integrable quantum systems [31,32,33,34,35]. A model of particular interest has been the S = 1/2 XXZ chain. At the SU(2) symmetric point (the XXX chain), the low T properties of this antiferromagnet are expected to be in the universality class of L in (1.1) at θ = π. So studies of the XXX model will explicate the nature of spin transport at temperatures T < T 0 at θ = π. This is a regime for which our paper has no results (although, spin transport in the low T < ∆ regime for θ = 0 was studied in Ref [4]). However, we have examined the higher temperature regime T 0 < T < T (2) max in Section III B, and found that our numerical results are not inconsistent with the presence of spin diffusion. It would then seem natural that diffusion may also exist for T < T 0 , although, this is, of course, not a rigorous argument. The latest numerical evidence for the XXX chain [33,35] seems to be consistent with the existence of diffusion.
Here ,we wish to issue a small caution towards the method used to study transport in Refs [31,32,33,35] (this caution does not apply to Ref [34]). These works computed a 'stiffness', which is the co-efficient of zero frequency delta function in the frequency dependent conductivity. By its very construction, such a quantity is defined at ω = 0 in a finite system; so implicitly, the limit ω → 0 has been taken before the L → ∞ limit. The considerations of this section make it amply clear that such a procedure is potentially dangerous.

V. IMPLICATIONS FOR EXPERIMENTS
First, we summarize the main theoretical results of this paper. We have shown that there is an intermediate temperature range over which the static and dynamic properties of a large class of one-dimensional Heisenberg antiferromagnets are described by the deterministic, continuum model defined by (1.14) and (1.18). For p-leg ladders of spin S ions, this temperature range rapidly becomes quite wide as Sp increases. For 2Sp even, our universal results hold for ∆ < T < T (2) max , where ∆ is the ground state energy gap, and T (2) max is estimated in (1.23), while for 2Sp odd, they hold for T 0 < T < T (2) max , where T 0 is an energy scale measuring the strength of logarithmic corrections at the lowest T . Both T 0 and ∆ become exponentially small as Sp increases, and so the intermediate temperature regime is clearly defined. The dynamical properties of such antiferromagnets are encapsulated in the scaling forms (1.21), which relate them to universal functions dependent only upon two thermodynamic parameters: the antiferromagnetic correlation length ξ(T ), and the uniform spin susceptibility χ u (T ). We obtained information on the universal functions in Sections II and III, including exact results on the spin-wave damping, while exact results for ξ(T ) and χ u (T ) were presented in Section I.
We now briefly review the work in the 70's on the dynamics of classical antiferromagnetic chain (Refs [9,10] and references therein). As we saw in Section I, there is a window of temperatures T (1) max < T < T (2) max (T (1) max was estimated in (1.3)), over which (1.24) and (1.25) are valid, where our results apply to purely classical models; so there is a common regime of validity between our and earlier work. These earlier classical results were all obtained on studies of lattice antiferromagnets, and all used some variant of the short-time moment expansion to extrapolate to the long time limit by a physically motivated ansatz e.g. the memory function formalism; however, there is a degree of arbitrariness in any such ansatz. In a regime where their correlation length ξ ≫ a (where a is a lattice spacing), and the wavevectors ka ≪ 1, their results should be described the continuum model we have discussed here. However, we expect their short-time methods to be exact only for kξ ≫ 1. Indeed, our paper is the first to provide a proper description of the scaling structure, along with quantitative information on the scaling functions, in the non-perturbative regime kξ ≪ 1. Consistent with these expectations, we saw in Section II that our result (2.20) for the spin-wave damping for kξ ≫ 1, and T (1) max < T < T (2) max , was in precise agreement with that of Reiter and Sjölander [10]. It should be noted, however, that our result (2.19) for the damping remains exact over a much wider window of temperatures (T 0 , ∆ < T < T (2) max ), including when χ u and ξ have the quantum renormalized T -dependence in (1.9) and (1.10). Moreover, we believe, despite conjectures by Reiter and Sjölander to the contrary, that the results in Ref [10] are not exact for kξ ≪ 1-we saw in Section III B that our value for the spin diffusion constant (assuming the existence of diffusion), D, disagreed with theirs.
Among static properties, neutron scattering experiments [38] have measured the correlation length, ξ(T ), and these have been compared to purely classical theories in which ξ(T ) behaves like (1.25). At lower temperatures, ξ(T ) should exhibit the logarithmic temperature dependence in (1.10), arising from quantum fluctuations. Combined with measurements of the uniform susceptibility, χ u (T ), a rather precise test of the quantum-renormalized static theory should then be possible.
Dynamic tests of the theory have focused mainly on the linewidth of the spin-wave excitations in the regime kξ ≫ 1. The measured linewidths have been compared [39,41] with the prediction of the classical theory [10], which yields the result (2.20). This is in general agreement with the theory, but a quantitative discrepancy was observed [41] for S = 3/2. We think it would be useful to compare the experiments with our new exact result (2.19), while using the actual experimentally observed values of ξ and χ u .
We think future neutron scattering experiments should also examine the interesting regime kξ ≪ 1. Here we have provided new, (numerically) exact results in Section III. In particular, there is some interesting physics in the structure of the frequency-dependent lineshapes in Figs 3 and 5, and these should be subjected to experimental tests. Also, we can easily generate additional universal spectra at other positions in the energy-momentum space, as needed.
The results in Section III also provide quantitative predictions for NMR experiments on spin chains. The nuclear relaxation rate 1/T 1 is given by local, low frequency, dynamic structure factor of the electronic spins. This has two contributions, one from the ferromagnetic component given by the correlator of L, and the other from the antiferromagnetic component given by the correlator of n. Let us parameterize the electronic spin S i by where a is the lattice spacing. If we assume that where Γ is related to the hyperfine coupling, and ω N → 0 is the nuclear Larmor frequency. The electron spin correlator has to be evaluated in the presence of an applied magnetic field H, and the electron Larmor precession can usually be neglected. However, for the case where there is spin diffusion, as in the assumed form (3.4), then this electron precession must be included for the Fourier transform is divergent at low frequencies. Combining (5.1) and (5.2) with the results of Section III, we obtain where D is estimated in (3.5,3.7). If we ignore logarithmic factors, the first antiferromagnetic term in (5.3) is of order Γ/T while the second ferromagnetic term is of order (Γ/T )(T a/c) 2 (T /H) 1/2 ; either term could be dominant, depending upon the magnitude of H. Further, our result (5.3) has assumed the existence of spin diffusion, but we expect that (5.3) will provide a reasonable quantitative estimate of the H and T dependence for experimental purposes, even if this assumption is not entirely correct in its details: there is clearly a long-time tail in Fig 7, even if it is not precisely diffusive. Finally, we compare our theoretical predictions with quantum Monte Carlo simulations on odd leg ladders. Numerical results for χ u (T ) have been obtained recently by Frischmuth et al. [42] on antiferromagnets with S = 1/2 and p = 3, 5. We compared their results for p = 5, with our new result (1.9,1.12): this is shown in Fig 12. The fitting parameters in this comparison are the values of c and T 0 . There is an arbitrariness in choosing the ranges of T over which to fit the intermediate T prediction (1.9), and this can lead to some variation in the values of c and T 0 . A reasonable set of values are c = 2.06Ja and T 0 = 0.0058J, are used in Fig 12. The value of c is roughly consistent with that estimated earlier [42], but the value of T 0 does appear to be rather small. This appendix will derive the relationship (1.12) between the two energy scales T 0 and Λ M S associated with the quantum O(3) non-linear sigma model at θ = π. This model has a flow from the high energy fixed point at g = 0 to the low-energy fixed point at g = g c . The flows near both fixed points are marginal: Λ M S is the energy scale characterizing the flow away from the the g = 0 fixed point, while T 0 is the scale characterizing the flow into the g = g c fixed point. These scales appear in the logarithmic corrections that appear in both the high temperature (Eqns (1.9,1.10)) and low temperature (Eqn (1.7)) limit. A complete Bethe ansatz analysis of the flow between the two fixed points matches the two scales, and leads to the relationship (1.12).
We will perform the matching by considering various limiting regimes of the free energy density, F Q , of Z Q as a function of T and H. This is will allow us to make an intricate series of mappings between numerous results which have appeared recently in the literature.
First, let us consider the low T and low H regime, where T, H ≪ T 0 . In this regime, the model is in the vicinity of the g = g c fixed point, which is the k = 1, SU(2) Wess-Zumino-Witten model. We know from (1.7), that for H ≪ T ≪ T 0 , the free energy has a contribution From this we anticipate that for T ≪ H ≪ T 0 , we will have a corresponding contribution where C 1 is a universal number we would like to determine. The universality of C 1 implies that we can use any model which is in the vicinity of the k = 1, SU(2) Wess-Zumino-Witten fixed point. In particular, Lukyanov [22] has recently computed the detailed H and T dependence of the free energy of the S = 1/2 antiferromagnetic chain with nearest neighbor exchange which fulfills this requirement, and we can use his results to obtain C 1 . In particular, from (3.18) and (3.20) of Ref [22] we determine that for the nearest neighbor S = 1/2 antiferromagnet, and that We emphasize that the result (A3) is non-universal, while (A4) is universal i.e. only the latter is a property of the continuum O(3) non-linear sigma model at θ = π. In the result (A4) we have set T = 0, as this is the limit in which we shall use in the following.
Let us now consider the vicinity of the g = 0 fixed point for non-zero H at zero temperature i.e. for H ≫ Λ M S ∼ T 0 at T = 0. In this case, a renormalized perturbation theory in g can be used to determine the free energy density, and this was carried out by Hasenfratz et al [21]. They obtained Finally, we need to match (A4) and (A5) by computing the T = 0 free energy for H = 0 in both the low and high field limits of the O(3) non-linear sigma model at θ = π. Fortunately, precisely this computation was carried out by Fateev, Onofri and Zamolodchikov [23]. They obtained (A5) in the limit H ≫ Λ M S , while the result for the opposite limit H ≪ Λ M S is in (4.97) of Ref. [23]: Comparing (A6) with (A4), we immediately obtain (1.12).

APPENDIX B: SHORT TIME EXPANSION: DETAILS
We start with the expression (2.2) and insert (2.4) into it and the expand to first order to find where Here we have dropped some additive constants and integrated by parts in places. As advertised before, the magnetic field in the z direction parametrised by m 2 adds a mass term to the action and makes the φφ * propagator infrared finite. The last term in the action arises from the Jacobian of the delta functionals δ(n · L)δ(n 2 − 1) in the measure. This Jacobian is infinite for a continuum system. To regularize it, we can introduce a discrete lattice in space. In that case the parameter a is the lattice constant, or equivalently the volume of the Brillouin zone.
Using (2.4) in the dynamical equations (2.1), we find the following equations of motion: These equations differ from (2.5) in that we have absorbed a factor of c = ξ 1/2 into the definition of time so that time and distance have the same units.
We have here written the nonlinearities as the inhomogeneous part of a set of first order equations. The solution to the initial value problem is given by The columns of the 2 × 2 matrix K(t) is made up of the two linearly independent solution vectors of the homogeneous problem.
Armed with this, we can solve for the fields in terms of the initial conditions and iterate the solution to go to higher orders by plugging the solution back in A and B. Also calculations simplify very much if we work out everything back in real space. So at last we write down the final iterative form from which the entire perturbation series can be generated.
To simplify notation let us denote the zeroth order solutions for the fields as ψ (0) (x, t) and φ (0) (x, t).
Now we can proceed to evaluate each of the correlation functions in (B4). As an example let us look at ψ(x, t)ψ * (0, 0) . The correlations to zeroth order are To compute the one loop correlation ψ(x, t)ψ * (0, 0) , we write down ψ(x, t)ψ * (0, 0) = ψ (0) (x, t)ψ (0) * (0, 0) On the second and third terms of the RHS, we can replace φ(ψ) by φ (0) (ψ (0) ) since it is already first order in 1/ξ. Now all the correlation functions on the RHS are expressed in terms of correlation functions at time t = 0. So we evaluate them using the partition function (B1). Other correlations can be evaluated in the same manner. So we shall only write down the final answers here.
Adding them all and taking the limit m → 0 leads to the expression quoted in (2.6). Note that all terms which diverge as m → 0, cancel each other only when we evaluate the O(3) invariant correlation L(x, t) · L(0, 0) .

APPENDIX C: EVALUATION OF THE ONE-PARTICLE DISTRIBUTION FUNCTION
Let us start with equation 4.8. Following Jepsen's notation we write with the definitions that The second term in the RHS of (C1) clearly represents the probability that the zeroth particle stays in the zeroth trajectory at time t, while the first term represents the probability that it has been scattered to another trajectory. The second term can be easily simplified using the definitions given above and in (4.5). It is (δ(y + ct) + δ(y − ct))P (1) Here p is defined so that 0 < 2ct − pL < L. As mentioned before, the peaks at times T , 2T , 3T , . . . are due to the exact recurrence of some configurations. We plot this function again for N = 50 and 51, L = 1, c = ±1 in Fig  11. Note that the peaks at times T , 2T . . . disappear. The first term will be handled in a manner similar to above. For simplicity we shall evaluate only the autocorrelation part. To do this, we set y = nL where n is an integer.
Using previous definitions we find that Here it will be convenient to define an integer p as before such that 0 < ct − pL < L.
Carrying out the integrals above we get The above function P (2) (t) has been plotted in Fig 10.