Solution of the classical Yang--Baxter equation with an exotic symmetry, and integrability of a multi-species boson tunneling model

Solutions of the classical Yang-Baxter equation provide a systematic method to construct integrable quantum systems in an algebraic manner. A Lie algebra can be associated with any solution of the classical Yang--Baxter equation, from which commuting transfer matrices may be constructed. This procedure is reviewed, specifically for solutions without skew-symmetry. A particular solution with an exotic symmetry is identified, which is not obtained as a limiting expansion of the usual Yang--Baxter equation. This solution facilitates the construction of commuting transfer matrices which will be used to establish the integrability of a multi-species boson tunneling model. The model generalises the well-known two-site Bose-Hubbard model, to which it reduces in the one-species limit. Due to the lack of an apparent reference state, application of the algebraic Bethe Ansatz to solve the model is prohibitive. Instead, the Bethe Ansatz solution is obtained by the use of operator identities and tensor product decompositions.


Introduction
The Yang-Baxter equation is a cornerstone in the theory of integrable quantum systems. It facilitates the construction of commuting transfer matrices which generate sets of conserved operators. A limiting procedure may be performed to obtain the classical Yang-Baxter equation [3]. Despite the name, this equation also has a long history in the study of integrable quantum systems [11,15,16]. In some cases these models have been formulated through the notion of a Gaudin algebra [2,12,21], which turns out to be an essentially equivalent approach up to generalisation of the classical Yang-Baxter equation [1]. Below, the use of a non-standard form of the classical Yang-Baxter equation will be employed, which affords solutions that do not have the property of being skew-symmetric. A particular case will be identified within a class of solutions of this equation, and analysed in detail.
The motivation for this study arises from a recent construction of an integrable generalisation, and Bethe Ansatz solution, of the p + ip-pairing Hamiltonian [19]. With the benefit of some hindsight, it will be shown below that this result can be re-derived in a more direct manner. Not only does this provide a clearer understanding of the integrability of the model derived in [19], it also opens up new avenues of applications. As an example, an integrable, multi-species, boson tunneling model will be derived. Following the experimental realisation of a two-species Bose-Einstein condensate in 2008 [35], there has been considerable effort to study such systems. Tunneling models for two-species systems were quickly formulated [8,22,25,34], and the general area of research continues to attract significant activity [6,7,20,37,39]. Having an integrable model for multi-species tunneling offers an opportunity for a different perspective on tunneling phenomena.
The integrable model which will be derived below is described in terms of 2L mutually commuting sets of canonical boson operators, representing 2L degrees of freedom for the system. These operators are labeled by two indices, j = 1, ..., L denoting the boson species, and X = A, B denoting the two potential wells between which tunneling occurs. The operators satisfy the relations Number operators for each species in each well are given by The Hamiltonian is given by The coupling parameter U characterising the inter-boson interactions is positive for repulsive interactions and negative for attractive interactions. The well bias, such as that due to the presence of an applied external field, is parametrised by µ. The parameters E j are the amplitudes associated with each species for tunneling between the wells. The signs of E j are not important, as they can be changed by a unitary transformation b j,X → −b j,X . There are no constraints imposed on any of the coupling parameters, they are all free, real variables. When L = 1 (1) reduces to the two-site Bose-Hubbard model, which endures as the subject of widespread investigation, e.g. [4,5,13,24,26,40]. Thus (1) represents an L-species generalisation. This Hamiltonian commutes with each of the operators N j , j = 1, .., L, indicating that the number of particles for each species is conserved. The main objective below is to construct an additional L conserved operators, establishing that the model is integrable.
In Sect. 2 a general Lie-algebraic approach will be taken to construct a family of commuting transfer matrices. It will be seen that the classical Yang-Baxter equation emerges in a natural way, in a form that does not assume skew-symmetry. This form has been promoted by Skrypnyk in a series of works [27,28,29,30,31,32,33]. Relationship to the usual form of the Yang-Baxter equation will be discussed. A study of some provisional elements of the associated representation theory will lead to the identification of a distinguished solution of the classical Yang-Baxter equation with an exotic symmetry, in a sense which will be defined. In Sect 3. attention turns to the construction of the model, and its solution. Due to the lack of an apparent reference state, application of the algebraic Bethe Ansatz to solve the model in general is prohibitive. Instead the Bethe Ansatz solution will first be derived in the limiting, but non-trivial, case when N j = 1 for all j = 1, ..., L. Based on this result the general solution will be obtained through a tensor product construction. Concluding remarks are offered in Sect. 4.

Commuting transfer matrices from the classical Yang-Baxter equation
Let {e j k : j, k = 1, ..., n} denote the standard basis elements for End(C n ). Introduce an abstract, infinite-dimensional complex Lie algebra G with generating functions and set For u, v ∈ C, let r(u, v) ∈ End(C n ⊗ C n ), which is expressible as The symmetry properties where t 2 denotes the partial transpose taken in the second component of the tensor product, are assumed. These assumptions provide for a different approach than those followed in [27,28,29,30,31,32,33]. The commutation relations imposed on G are compactly given by where the subscripts denote the actions on tensor components of C n ⊗ C n . Assuming that the Jacobi identity holds, such that an associative enveloping algebra can be constructed in the usual way, it follows from the commutation relations (4) that the transfer matrix This result does not depend on any properties of r(u, v) other than (2,3), as well as (4) being compatible with the Jacobi identity. It can be shown that a sufficient condition for the Jacobi identity to hold is the imposition that the non-standard classical Yang-Baxter equation (CYBE) holds. Recall that for the complex vector space C m , a Lie algebra representation π : G → End(C m ) is required to preserve the commutator in the sense The imposition of the CYBE also provides the notion of an adjoint representation π ad : G → End(C n ). The identification satisfies (7) as a consequence of (4) and (6). The above construction is non-standard in the sense that there is no assumption that the solution of (6) possesses the skew-symmetry property Consequently such solutions lie outside the Belavin-Drinfel'd classification of solutions [3]. The skew-symmetry property arises naturally in cases where a solution of (6) is obtained from a solution R(u, v) ∈ End(C n ⊗ C n ) of the Yang-Baxter equation with the unitarity condition where I is the identity operator. If R(u, v) depends on some parameter η, and admits an expansion (10) while it follows from (9) that It is then seen that (10) and (11) combined are equivalent to (6).
As an example, let P : C n ⊗ C n → C n ⊗ C n denote the permutation operator with the action A simple solution of (6) is This solution has the property (11), so it also provides a solution of (10). Hereafter set n = 2. One family of solutions of (6) which does not generally possess the skew-symmetry property (11), but does satisfy (2,3), was presented by Skrypnyk [30]. This solution can be represented in matrix form as where c ∈ C is arbitrary. The above solution is skew-symmetric only for the particular choice c = 0. In this instance the algebra G coincides with the trigonometric/hyperbolic sl(2) Gaudin algebra [21,36]. The primary interest here is in the case c = 1 which, as will be explained below, gives rise to an exotic symmetry. For this choice This solution may also be expressed as . This form is convenient for calculations that follow in Subsect. 2.1.

Symmetries
Symmetries of solutions of (6) may be characterised in terms of the onedimensional representations π 0 of their associated Lie algebra G, i.e. π 0 : G → C. From (7) it is seen that for any one-dimensional representation Equivalently, a one-dimensional representation may be expressed as a mapping ρ(T (u)) = J (u) ∈ End(C 2 ) satisfying In [32] such a solution is termed a generalised shift element.
For the case of the solution (12), it is seen that all constant elements of End(C 2 ) satisfy (15). It is therefore natural to identify the symmetry of (12) as gl (2). On the other hand for generic values of c in (13), solutions of (15) are given by constant diagonal elements of End(C 2 ). In this instance the symmetry of (13) is identified as u(1). For the specific choice c = 1 in (13) leading to (14), an exotic symmetry emerges. Setting J (u) = D + uO where D is diagonal and O is off-diagonal it is found that (15) is satisfied by J (u). First where simplification to the final expression has relied on the fact that σO = −Oσ. Now every diagonal operator on C 2 may be expressed as D = aI + bσ for some a, b ∈ C. It then follows that Because this solution for J (u) has an explicit dependence on the parameter u, it is considered to be exotic compared to the gl(2) symmetry of (12) and the u(1) symmetry of (13). For later use the notation G will be adopted for the Lie algebra associated with the case of the solution (14).

General representations
General representations of G are first obtained by an evaluation homomorphism ρ z to the enveloping algebra U(sl (2)) of sl (2). Adopting the convention for the sl(2) generators {S z , S + , S − } to satisfy the commutation relations preserving the commutation relations (4) for the solution (14). Every representation of sl(2) then affords a a representation of G. For the choice of the spin-1/2 representation of sl(2) this coincides with the adjoint representation (8).
A second evaluation homomorphism to the enveloping algebra U(sl(2)), denoted ρ (0) , is deduced from the one-dimensional representation. Specifically, where α, β, γ, δ ∈ C are arbitrary, which again preserves the commutation relations (4) for the solution (14). Taking the tensor product of these evaluation homomorphisms provides evaluation homomorphisms acting on tensor product of copies of U(sl (2)). The general map It then follows that τ (u) = ρ ⊗L (t(u)) ∈ U(sl(2)) ⊗L satisfying is a transfer matrix expressed in terms of spin operators.
Hereafter set δ = −α, and γ = β with α, β ∈ R, giving where C = 2(S z ) 2 + S + S − + S − S + is the sl(2) Casimir element and With the restriction z j ∈ R, j = 1, ..., L, which is imposed from now on, the {T j : j = 1, ..., L} form a set of mutually commuting, self-adjoint operators. Consequently they are simultaneously diagonalisable. It is seen that the union of this set, the identity operator, and the set of Casimir elements {C j : j = 1, ..., L} is the full set of functionally independent operators that can be extracted from the transfer matrix. Matrix forms of transfer matrix are obtained by choosing a representation of sl(2) for each copy in the tensor product U(sl(2)) ⊗L . It should be noted that these representations need not be isomorphic, and throughout the parameters α, β, z 1 , z 2 , ..., z L are arbitrary real parameters. This provides a range of potential applications for constructing integrable models.

Integrable, multi-species, boson tunneling model
Having identified the exotic symmetry of the solution (14), and outlined the construction of representations of G, attention turns towards an application through the construction of an integrable, multi-species, boson tunneling model. The Hamiltonian is defined as Next the spin operators are represented in terms of boson operators through the Jordan-Schwinger map where the spin s j of each representation is determined by 2s j = N j . This maps (17) to (1) with the identification µ = 2Uα and E j = 2Uβz j . For simplicity it will be assumed hereafter that the z j , and consequently the E j , are distinct parameters.
Since the operators (16) commute with (17), application of the Jordan-Schwinger map to the T j yields a set of L operators which necessarily commute with (1). Along with the N j , which also commute with (1), the construction for the model is one which establishes that it is integrable since the number of independent conserved operators is the same as the number of degrees of freedom. Closely related to the property of integrability is the notion of an exact solution. The next step is to derive an explicit form for such a solution.

Bethe Ansatz solution
Below it will be shown that the eigenvalues Λ j of the conserved operators (16) have the form such that the parameters {v n : l = 1, ..., N}, known as the Bethe roots, satisfy the Bethe Ansatz equations Consequently, the energy eigenvalues E for the Hamiltonian (1) have the form Let the representation space associated with highest weight s be denoted V (s), such that for a system with N j bosons of each species j the Hilbert space is The tensor product state consisting of all the highest weight states is not an eigenstate of the transfer matrix, except in the case when β = 0. Consequently, the algebraic Bethe Ansatz method of solution for calculating eigenvalue expressions of the operators (20) cannot be easily applied due to the absence of an apparent reference state. Here a different approach is taken, based on operator identities in the first instance for a specific case. Using the method of proof by induction, the result will be extended to the general case by examining tensor product decompositions.

N-species case
First the case where N j = 1 for all j = 1, ..., L, so N = L, will be considered. The quasi-classical limit of the Off-diagonal Bethe Ansatz approach [38] was utilised in [19] to obtain the solution in this instance. While the physical interpretation of the p+ip pairing Hamiltonian adopted in [19] is a completely different context to the one considered here, at the mathematical level it is entirely equivalent and can be seen as the case where all copies of sl (2) algebras are represented by spin 1/2. It was found that the eigenvalues λ j are given by (18) subject to the Bethe Ansatz equations (19) with the specialisation s j = 1/2 for all j = 1, ..., L. In [9] a method which generalises the algebraic Bethe Ansatz was used to confirm the result.
The approaches of [9,19] are quite technical. Below a different approach will be taken, following [17], which reproduces the equations in a streamlined manner. To simplify calculations, first define the operators which, in the spin-1/2 case, only differ from (16) by constant terms. The starting point uses a result given in [9] that the operators T j satisfy the quadratic identities Since these are operator identities it follows that the eigenvalues λ j of the operators T j satisfy analogous relations: Define Q(u) to be a polynomial of order L satisfying Since this is a homogeneous linear system of L equations for the L + 1 coefficients of the polynomial Q(u), there exists a non-trivial solution. Set which assumes that the roots of Q(u) do not have multiplicities. The validity of this assumption will be reviewed later. Now, provided Q(z 2 j ) = 0 for all j = 1, ..., L, then On the other hand from (21) where For equality of the expressions (22) and (23) it is required that Set which is a polynomial of order L. It follows from (24) that S(z 2 j ) = 0, j = 1, ..., L which, along with the consideration of the asymptotic behaviour of (25) as u → ∞, establishes that Evaluating S(v n ) through (25) and (26) and equating these expressions then yields the Bethe Ansatz equations which are equivalent to (19). It can be verified that (18) holds.
It needs addressing at this point that the assumptions that the roots of Q(u) do not have multiplicities, and that Q(z 2 j ) = 0, are not always valid. In the β = 0 limit numerical solutions of the Bethe Ansatz equations confirm that both of these situations occur [10,14,18,23]. Two Bethe roots may coalesce at a point z 2 j through which they transition between being realvalued and a complex-conjugate pair. Moreover in this limit some of the L Bethe roots must be seen as being infinite-valued, and there are well-identified points at which several additional roots collapse to be zero-valued or infinitevalued. All of these instances may be interpreted as singular solutions of (19), but can be simply regularised by expressing (27) in a polynomial form So, while the above derivation of the Bethe Ansatz equations breaks down when the roots of Q(u) have multiplicities, or Q(z 2 j ) = 0, indicators from available numerical studies [10,14,18,23] suggest that the equations are nonetheless valid.

The general case
The validity of the formulae (18,19) in general follows by an inductive argument on each of the species number labels N j , starting with the case N j = 1 for all j = 1, ..., N which has been established above. Every finitedimensional irreducible sl(2) module is contained in some tensor product of spin-1/2 modules, and in particular where N = N 1 + ... + N L as before. Now consider a system of M ≥ L species, for which it is determined that for j = M − 1, M lim z M →z M −1 T j =2(S z j ) 2 + 2αS z j + βz j (S + j + S − j ) whereby the tilded notation refers to a new system of N particles of M − 1