Exact physical quantities of a competing spin chain in the thermodynamic limit

We study the exact physical quantities of a competing spin chain which contains many interesting and meaningful couplings including the nearest neighbor, next nearest neighbor, chiral three spins, Dzyloshinsky-Moriya interactions and unparallel boundary magnetic fields in the thermodynamic limit. We obtain the density of zero roots, surface energies and elementary excitations in different regimes of model parameters. Due to the competition of various interactions, the surface energy and excited spectrum show many different pictures from those of the Heisenberg spin chain.


Introduction
Quantum integrable models [1] are very important to analyze some non-pertubative properties of quantum field/string theory [2,3]. Moreover, the exact solutions and physical properties of these models can provide the strict benchmarks for many important physics issues, and sometimes it can exactly predict and explain the results of experiments [4][5][6]. In recent years, the study of quantum integrable models play an important role in the non-equilibrium statistical physics [7][8][9][10], condensed matter physics [11], cold atom physic [12,13], superstring theory AdS/CFT [14][15][16] and so on.
For the integrable models with U(1) symmetry, the exact solutions of the models can be obtained by the conventional Bethe ansatz. In addition, due to the homogeneous Bethe ansatz equations (BAEs) and the regular pattern of the Bethe roots, the thermodynamic properties can be directly calculated by the thermodynamic Bethe ansatz (TBA) [17,18].
When the U(1) symmetry of integrable systems is broken, the off-diagonal Bethe ansatz can be used to solve the systems based on the algebraic analysis [19]. However, since the exact solutions of the systems are described by the inhomogeneous T − Q relations [20,21] and the resulting inhomogeneous BAEs have the inhomogeneous term, the pattern of Bethe roots is not clear and the TBA method can not be applied. Recently, a novel Bethe ansatz scheme has been proposed to calculate the physical quantities of quantum integrable systems with or without U(1) symmetry [22,23]. The key point of the scheme lies in parameterizing the eigenvalue of transfer matrix by its zero roots instead of the Bethe roots. Through this method, the homogeneous BAEs and the well-defined patterns of zero roots can be obtained. Based on them, the thermodynamic properties and exact physical quantities of the systems in the thermodynamic limit can also be calculated. In this paper, we study an isotropic quantum spin chain which includes the nearest neighbor (NN) [24], next nearest neighbor (NNN) [25], Dzyloshinsky-Moriya (DM) interactions [26,27], chirality three-spin couplings [28] and unparallel boundary magnetic fields [29]. The density of zero roots, surface energy and elementary excitations in different regimes of model parameters are obtained.
The paper is organized as follows. Section 2 serves as an introduction to the model and explain its integrability. In section 3, we give the patterns of zero roots in the different regimes of model parameters. In section 4, we calculate the surface energies induced by the boundary magnetic fields. In section 5, we study the typical bulk elementary excitations in the system. The boundary excitations are computed in section 6. In section 7, we calculate the surface energies in ferromagnetic regime. Concluding remarks are given in section 8. A simple method is introduced in Appendix A..

Integrability of the model
The model Hamiltonian reads Here H bulk describe the interactions in the bulk which includes the NN, NNN and chiral three spin couplings with the form of where σ α j (α = x, y, z) is the Pauli matrix along the α-direction on the j-th site, and 2N is the number of sites. We note that the convention σ 2N +1 = 0 has been used. H L quantifies the left boundary terms which includes the boundary magnetic field along the z-direction and the anisotropic and DM interactions of the first bond where p is the strength of magnetic field, a 2 and ap quantify the spin-exchanging and DM interactions respectively, and D z 1 is the unit vector along the z-direction. H R characterizes the right boundary terms which includes the boundary magnetic field lies in the x − z plane, anisotropic and DM interactions of the last bond also constrained in the x − z plane. Thus where q and ξ are the boundary parameters, D x 2N is the unit vector along the x-direction and D z 2N is the unit vector along the z-direction. We should note that the boundary fields are unparallel boundary and the U(1) symmetry of the system are broken. The hermitian of the Hamiltonian (2.1) requires that the model parameter a is pure imaginary and the boundary parameters p, q, ξ are real. Moreover, the integrability of the system (2.1) requires that the couplings J 1 , J 2 , J 3 satisfy the relationships where the index j is the summation index in H bulk (2.2). The Hamiltonian (2.1) is constructed by using the R-matrix and the reflection matrices K ± based on the quantum inverse scattering method. The R-matrix defined in the tensor space where u is the spectral parameter and P 1,2 is the permutation operator. The R-matrix (2.7) satisfies the quantum Yang-Baxter equation (QYBE), which satisfies the reflection equation (RE) where R 2,1 (u) = P 1,2 R 1,2 (u)P 1,2 . The dual reflection matrix K + 1 (u) is satisfying the dual reflection equation The monodromy matrix T 0 (u) and the reflecting oneT 0 (u) are constructed by the R-matrices where V 0 is the auxiliary space, ⊗ 2N j=1 V j is the quantum space, and {θ j |j = 1, · · · , 2N} are the inhomogeneity parameters. The transfer matrix t(u) is defined as where tr 0 means the partial trace over the auxiliary space. The Hamiltonian (2.1) is generated by the transfer matrix as where The QYBE (2.8), the RE (2.10) and its dual (2.12) guarantee the integrability of the model described by the Hamiltonian given by (2.1). Moreover, using the properties of the R-matrix one may easily prove that t(u) = t(−u−1) and the following operator identities [19] t(θ j + a)t(θ j + a − 1) = a(θ j + a)d(θ j + a − 1), j = 1, · · · , 2N, From the definition (2.14), we know that the transfer matrix t(u) is a polynomial operator of u with the degree 4N + 2. Denote the eigenvalue of the transfer matrix t(u) as Λ(u). From above analysis, we know that the eigenvalue Λ(u) satisfies Λ(θ j + a)Λ(θ j + a − 1) = a(θ j + a)d(θ j + a − 1), j = 1, · · · , 2N. (2.22) Obviously, Λ(u) is a degree 4N + 2 polynomial of u and can be parameterized as where {z j |j = 1, · · · , 2N +1} are the zero roots of the polynomial. Putting the parameterizing (2.23) into (2.22), we obtain the BAEs The above 2N equations and (2.21) can determine the 2N + 1 unknowns {z j } completely.
(2.25) and (2.26) can determine the 2N +1 zeros roots {z j } in the homogeneous limit in finite system size. Moreover, the energy spectrum of the Hamiltonian (2.1) can be determined by the zero roots as where the function a n (u) is given by a n (u) = 1 2π

Patterns of zero roots
We first study the solutions of zero roots {z j } at the ground state. For convenient, we choose all the inhomogeneity parameters to be imaginary, {θ j ≡ iθ j }, and let {z j ≡ −iz j }.
In addition, we set the boundary parameters as p > 0 andq = q(1 + ξ 2 ) − 1 2 . From the numerical calculation and algebraic analysis, we find that the distribution of thez-roots at the ground state can be divided into following six different regimes in the upper p −q plane, as shown in Fig.1.   The blue asterisks indicatez-roots for {θ j = 0|j = 1, · · · , 2N} and the red circles specifȳ z-roots with the inhomogeneity parameters {θ j = 0.1(j − N − 0.5)|j = 1, · · · , 2N}.
We also find that the choice of pure imaginary inhomogeneities {θ j } does not change the patterns of zero roots {z j } but the roots density, as shown in

Surface energy
Now, we consider the surface energy induced by the boundaries. The surface energy is where E g is the ground state energy of present system and E p is the ground state energy of the corresponding periodic chain. In the thermodynamic limit, the distribution ofz-roots can be characterized by the density ρ(z). Furthermore, we assume that the density of inhomogeneity parameters 1/[2N(θ j −θ j−1 )] has the continuum limit σ(θ).
In regime I, substituting the corresponding pattern ofz-roots into BAEs (2.24) and taking the logarithm of the absolute value, we have ln |4|+ whereā = −ia. In the thermodynamic limit, we assume that the zero roots and inhomogeneities have continuum densities .
Taking the continuum limit of Eq. (4.1) and replacingθ j with λ, we obtain where b n (λ) = 1 2π 2λ λ 2 +n 2 /4 . Eq.(4.2) is a convolution equation and can be solved by the Fourier transformation. The solution ofz-roots density is whereb n (k) = sign(k)ie −|nk| . From now on, we use σ(θ) = δ(θ). In the thermodynamic limit, α tends to infinity. The ground state energy of the Hamiltonian (2.1) can thus be expressed as whereã n (k) = e −|nk| is the Fourier transformation of a n (λ). The ground state energy of the system with periodic boundary condition can be obtained similarly. After tedious calculation, we obtain the surface energy in the regime I as dk. In the regime II, taking the logarithm then the derivative of the absolute value of BAE The Fourier transform gives (4.10) Then we obtain the surface energy in this regime as where e b (p), e b (q) and e b0 are given by Eqs.  Fig.4(a), we see that the surface energy of Heisenberg spin chain is smaller than zero, and is monotonically increasing with the increasing of |p|.
When p = 0, the surface energy is divergent, this is because that the strength of boundary magnetic field is quantified by 1/p. The results are similar to the those of the Heisenberg spin chain with parallel boundary fields [30,31]. While for the present model with a = 0, the surface energies can be larger or smaller than zero, and have two peaks and three minimums at some special values of |p|. At the point of p = 0, the surface energy arrives at its minimum.
The surface energy is smaller than that of Heisenberg spin chain if |p| is large, and is larger than that of Heisenberg spin chain if |p| is smalle.
The surface energies e b (p) with fixed a versus p are shown in Fig.4(b). Comparing

Figs.4(a) and (b)
, we find that if |p| is large which means that the boundary field is small, due to the existence of NNN, chiral three spin and DM interactions, the surface energy is smaller than that of the Heisenberg spin chain. We should note that the relation between e b (q) andq is the same as that between e b (p) and p, whereq = q/ 1 + ξ 2 .
The strength of boundary magnetic field along the z-direction is quantified by p or q up to a normalized scalar factor. The further numerical calculation of the analytical expression of surface energy shows that the curves of E b versus q are similar with those of E b versus p.
Thus we omit the figure of E b with the changing of q here. In Fig.4(c), we show the surface energies E b with given a versus the boundary parameter ξ. The ξ quantifies the twisted angle between two unparallel boundary magnetic fields, and quantifies the strength of magnetic field on the right boundary. If ξ is large, the twisted angle is large. At the same time, the right boundary magnetic field is small. From the blue dotted lines in Fig.4(c), which corresponds to the Heisenberg spin chain, we see clearly that if ξ is small, the magnetic field is strong thus the induced surface energy is large, as it should be. For the present system with a = 0, if ξ is small, the contributions of NNN, chiral three spin and DM interactions are large, which leads to the surface energy becomes small. Thus the behaviors of surface energies with a = 0 and a = 0 are totally different.
The surface energies e b0 versus the different values of parameter a are shown in Fig.4(d).
We note that the value of e b0 at the point of a = 0 is the surface energy of the Heisenberg spin chain with free open boundaries.
From above explanations, we conclude that the surface energy of present system is quite different from that of the Heisenberg spin chain.

Bulk elementary excitations
Next, we study the elementary excitations in the system. We first consider the excitations in the bulk. The bulk excitations in different regimes of boundary parameters are the same.
From the patterns of zero roots in the low-lying excited states, we find that the excitations can be characterized by breaking several conjugate pairs and putting the corresponding zero roots into the real axis, or the zero roots forming the conjugate pairs on the imaginary axis with more larger imaginary parts ± ni 2 (n > 2). Thus the system has two kinds of bulk elementary excitations. The first one is quantified by four finite real roots {±z 1 , ±z 2 } and the second one is quantified by two conjugate pairs {z n ± ni 2 , −z n ± ni 2 }, where the distribution of rest zero roots almost does not change and the related difference between ground and excited states can be erased by the rearrangement of Fermi sea in the thermodynamic limit.
As an example, we give the pattern of zero roots at the ground state (blue asterisks) and that at the first kind of excited sate (red circles) in the regime V with 2N = 8, which is shown in Fig.5(a). It is clear that there are four new real roots at the excited state. In the thermodynamic limit, the density difference δρ e 1 (k) between the ground state and the excited state is wherez 1 andz 2 can take arbitrary continuous values in the real axis. Thus the energy carried by this kind of excitation is δ e = δ e 1 (z 1 ) + δ e 1 (z 2 ), which covers the previous results obtained by using the conventional Bethe ansatz method for the periodic staggered (a = 0) spin chain [32]. The excited energies δ e 1 with given values of model parameter a versusz 1 are shown in Fig.5(b). From it, we see that the excited energy of the Heisenberg spin chain (a = 0) only has one peak at the point ofz = 0, while for the present model (a = 0), the excited energies have two peaks at finite ±z.
Now, we focus on the second kind of elementary excitation. In order to see the high strings (n > 2) excitations more clearly, we show the pattern of zero roots at the n = 3 excited state in Fig.6, where the ground state is still in the regime V. In the thermodynamic limit, the density difference δρ en (k) between the ground state and the excited state is wherez n is free. The related elementary excitation energy is (1 − e −|k| ) cosh(ak) (e −|(n+1)k|/2 + e −|(n−1)k|/2 ) cos(z n k) e −|k|/2 cosh (k/2) dk +2π(a n+1 (z n + ia) + a n+1 (z n − ia) − a n−1 (z n + ia) − a n−1 (z n − ia))] = 0, (5.4) which indicates that the bare contributions of the conjugate pairs with n > 2 to the energy is exactly canceled by that of the back flow of the continuum root density. Thus the conjugate pairs with n > 2 contribute nothing to the energy. However, the conjugate pairs do affect the scattering matrix among the real roots [33].

Boundary elementary excitations
Next, we consider the boundary excitations. Comparing with the zero roots distributions at the ground state, we find that the boundary excitations can exist in the regimes I-IV, where the boundary parameter − 1 2 < p < 1 2 or − 1 2 <q < 1 2 . The typical boundary excitation is  putting the boundary string from i(|p| + 1 2 ) to i( 1 2 − |p|), or from i(|q| + 1 2 ) to i( 1 2 − |q|). These two new boundary strings indeed are the solutions of BAEs (2.24) and would appear at the low-lying excited states.
As an example, we show the pattern of zero roots at the ground state (blue asterisks) and that at the excited state (red circles) with boundary string i( 1 2 − |p|) in the regime III with 2N = 8, which is shown in Fig.7(a). We can find in the excitation, the 4 roots at ±α and ±β of the ground state jump into the bulk string parts at ±i axes. The change of the zero roots ±α and ±iβ contribute nothing to the energy. Therefore, we omit the zero roots ±α and ±iβ in the following. The resulted density change δρ(k) between ground and excited states reads The corresponding excited energy is = −π(4a 2 − 1) · csc(π(|p| + a)) + csc(π(|p| − a)) . (6. 2) The excited energies δ ep with fixed values of a versus p are shown in Fig.7(b). From it, we see that the excited energy of present model is increasing with the increasing of boundary parameter |p| and has a minimum at the point of p = 0, which is very different from that of the Heisenberg spin chain. For the latter, the excited energy is decreasing with the increase of |p|.
We have computed the boundary excitations in other regimes and found that the excited energies has an unified form (6.2), although the resulted values are different. Please note that when considering the boundary excitations in the regime of − 1 2 <q < 1 2 , the p in Eq.(6.2) should be replaced by theq.

Surface energy in ferromagnetic regime
Furthermore, we study the surface energy in ferromagnetic regime. The corresponding Hamiltonian H f err is the negative of Hamiltonian (2.1), namely In region III (p ≥ 1 2 , 0 ≤q < 1 2 or 0 ≤ p < 1 2 ,q ≥ 1 2 ), all the zeros {z j |j = 1, · · · , N} are real as shown in Fig.8(a). Taking the logarithm then the derivative of the absolute value of BAE (2.24), we have The Fourier transform gives The ground state energy of the Hamiltonian (7.1) can thus be expressed as where the surface energy E f err b in this regime as  Fig. 8(b).

Conclusions
In this paper, we have studied the exact physical quantities of a competing spin chain including the NN, NNN, chiral three-spin couplings, DM interactions and unparallel boundary magnetic fields in the thermodynamic limit. We obtained the density of zero roots, surface can be presented in a simplified, explicit form, by doing the Fourier integral resulting in: δ e 1 (z) = −(4a 2 − 1) · π cosh(z + ia) + π cosh(z − ia) . This is solved uniquely for a given set of zeros z m in the physical strip by tanh resp. tan function (for any distribution of inhomogeneity parameters θ j ). Let us assume there are only two such zeros z 1 and z 2 , then l(u) = tan π 2 (u − z 1 ) + 1 2 π 2 (u − z 2 ) + 1 2 . (A.12) The shift + 1 2 is due to the convention (2.23). The logarithmic derivative and then inserting u = ±a and z m = iz m gives directly (5.2).
However, the method requires that there do not exist the zeros between the lines Re(z) = 0 and Re(z) = −1 at the ground state. For example, we can know that zeros of the ground state in ferromagnetic regime are mainly located in line Re(z) = − 1 2 4 from Section 7 . This will lead to an error in the Fourier transform (A.7).