Assembling integrable sigma-models as affine Gaudin models

We explain how to obtain new classical integrable field theories by assembling two affine Gaudin models into a single one. We show that the resulting affine Gaudin model depends on a parameter $\gamma$ in such a way that the limit $\gamma \to 0$ corresponds to the decoupling limit. Simple conditions ensuring Lorentz invariance are also presented. A first application of this method for $\sigma$-models leads to the action announced in [Phys. Rev. Lett. 122 (2019) 041601] and which couples an arbitrary number $N$ of principal chiral model fields on the same Lie group, each with a Wess-Zumino term. The affine Gaudin model descriptions of various integrable $\sigma$-models that can be used as elementary building blocks in the assembling construction are then given. This is in particular used in a second application of the method which consists in assembling $N-1$ copies of the principal chiral model each with a Wess-Zumino term and one homogeneous Yang-Baxter deformation of the principal chiral model.


Introduction
Constructing classical field theories which are integrable is difficult. Indeed, the property of integrability is established at the Lagrangian level if one is able to find a Lax connection whose zero curvature equation is equivalent to the equations of motion derived from the given action. This ensures, in particular, the existence of infinitely many integrals of motion in the theory. Yet finding such a Lax connection is somewhat of an art. Furthermore, in order to establish integrability at the Hamiltonian level one also needs to prove the existence of an infinite number of integrals of motion in involution. This is guaranteed if the Poisson bracket of the Lax matrix can be brought to a certain form. But doing so is also, in and of itself, an arduous task.
Faced with the difficulty of finding new integrable field theories from scratch, one approach to exploring the landscape of all possible Hamiltonian integrable field theories is to build new ones from old ones. One way of doing so is by deforming. For a very broad class of integrable field theories, the Poisson bracket of the Lax matrix is controlled by a certain rational function ϕ(z) of the spectral parameter z, called the twist function. As shown in [2,3], this twist function can be used at the Hamiltonian level to deform the theory in question, all the while preserving integrability. We may then perform an inverse Legendre transform to construct the actions of the resulting deformed field theories, which by construction will automatically be integrable.
It was recently shown by one of us in [4] that classical integrable field theories admitting a twist function can be seen as realisations of classical Gaudin models associated with untwisted affine Kac-Moody algebras. In fact, all of the examples considered in [4], which includes many integrable σ-models and their deformations, essentially correspond to affine Gaudin models with a single site. However, since the number of sites in an affine Gaudin model can be arbitrary, this opens up another avenue for constructing new integrable field theories from existing ones.
Indeed, the procedure detailed in the present article consists, roughly speaking, in taking existing integrable field theories which admit an affine Gaudin model description and 'placing' them at the individual sites of an N -site affine Gaudin model. This results in a new integrable field theory which couples these individual theories together. In fact, the strength of the coupling will be determined by the relative location of each site in the complex plane. Finally, since affine Gaudin models are intrinsically formulated at the Hamiltonian level, to complete the procedure one has to perform an inverse Legendre transform in order to construct the action for this new model. In [1] we presented the result of applying such a procedure to couple together N principal chiral fields on the same Lie group, each with a Wess-Zumino term.
To motivate why building an N -site affine Gaudin model in this way corresponds to coupling together the individual field theories, let us consider here the simpler situation of a Gaudin model associated with a finite-dimensional complex semisimple Lie algebra g [5]. Let u i for i = 1, . . . , N be a collection of N distinct complex numbers. Let I a(i) for a = 1, . . . , dim g denote the copy of a basis of the Lie algebra g that we think of as being formally 'attached' to the point u i , for i = 1, . . . , N . One may regard the basis elements I a(i) , i = 1, . . . , N of the N copies of g as the degrees of freedom of the model. They satisfy the Kostant-Kirillov Poisson bracket We may then express the quadratic Gaudin Hamiltonians (1.2) as H i = 1 2 res u i κ(L(z), L(z))dz. (

1.4)
A finite-dimensional integrable system is said to be a realisation of this Gaudin model if there is a morphism from the algebra of observables of the Gaudin model to that of the given integrable system. Moreover, the Hamiltonian of the integrable system should coincide, in the simplest case, with the image of some linear combination of the Hamiltonians (1.2) under this morphism. Concretely, the existence of such a morphism means that the integrable system should admit a Lax matrix of the same form as (1.3) and satisfying the same Poisson bracket relations as (1.3) does. In particular, the g-valued observables J (i) and J (j) should commute for i = j. Suppose now that we break up the collection of N sites {u i } N i=1 into two subsets of sizes N 1 and N 2 , respectively, which without loss of generality we take to be {u i } N 1 i=1 and {u j+N 1 } N 2 j=1 . Fixing a non-zero complex number γ = 0, we then choose to parametrise these two subsets of points as u i = z i − 1 2 γ −1 , i = 1, . . . , N 1 and u j+N 1 =ẑ j + 1 2 γ −1 , j = 1, . . . , N 2 . For convenience we also use the notationÎ a(j) := I a(j+N 1 ) andÎ (j) a := I (j+N 1 ) a for j = 1, . . . , N 2 . The collection of Hamiltonians (1.2) then naturally split into two subsets as well, namely for every i = 1, . . . , N 1 and j = 1, . . . , N 2 . We now observe that the first sums in each of the above expressions (1.5) is just the usual expression for the Gaudin Hamiltonians, as in (1.2), but corresponding to one of the two subsets of points. Moreover, the second sums on the right hand sides of (1.5) describe interaction terms between the degrees of freedom attached to the first N 1 sites and those attached to the last N 2 sites. It is clear from (1.5) that one can decouple these two sets of degrees of freedom by letting the 'coupling' γ go to zero. This corresponds to moving apart the first N 1 sites from the last N 2 sites while keeping all the relative positions within these two sets of sites fixed. In this way one can view the Gaudin model with N = N 1 + N 2 sites as the coupling of two separate Gaudin models, one with N 1 sites and the other with N 2 sites. It is important to note that all three Gaudin models are associated with the same Lie algebra g. In fact, by reversing the logic of the previous paragraphs, we see that the coupling of finitedimensional integrable systems which are realisations of Gaudin models associated with the same Lie algebra g is easy to perform at the level of the Lax matrix. Starting with the individual Lax matrices for both models, namely where J (i) 1 = I a ⊗ I a(i) for i = 1, . . . , N 1 and J (i) 2 = I a ⊗Î a(i) for i = 1, . . . , N 2 , the Lax matrix of the coupled system in (1.3) is obtained simply by adding the Lax matrices (1.6) together after suitably shifting their spectral parameters by some amount of the inverse γ −1 of the coupling parameter γ. Explicitly, we have L(z) = L 1 z + 1 2 γ −1 + L 2 z − 1 2 γ −1 . (1.7) The Hamiltonian of the coupled system, which includes all the interaction terms as in (1.5), is then constructed from this Lax matrix in the usual way as a linear combination of (1.4).
We now return to the field theory setting. Let g be the untwisted affine Kac-Moody algebra corresponding to g. It is useful to start by briefly recalling how classical integrable field theories which admit a twist function can be viewed as Gaudin models associated with g, referring the reader to [4] for the details.
The Gaudin model associated with the infinite dimensional Lie algebra g is defined in exactly the same way as in the finite-dimensional case recalled above, replacing g everywhere by g. The only subtlety is that, in doing so, the implicit finite sums over a = 1, . . . , dim g are now replaced by infinite sums over an index a labelling a basis of g, which one has to make sense of by working in a suitable completion of the algebra of observables. See [4] for a further discussion of this technical point which we shall not consider any further here. In particular, the Lax matrix of the affine Gaudin model takes on the exact same form as in (1.3) but where now the g-valued observables The former is precisely the twist function of the affine Gaudin model and we refer to the latter as the Gaudin Lax matrix. Among the many possible generalisations of the Gaudin model, one of particular interest for the present article involves including higher order poles at each u i , i = 1, . . . , N in the Lax matrix (1.3), or also higher order singularities at infinity by introducing polynomial terms. In fact, all the field theories we shall consider correspond to Gaudin models with the Lax matrices having at most double poles at each u i , i = 1, . . . , N and at infinity. Such Lax matrices can then also be rewritten in the form ϕ(z)∂ x + Γ(z, x) to read off the twist function and Gaudin Lax matrix.
So far the presentation we gave of Gaudin models in the affine setting has been completely parallel to the one given earlier in the finite case. However, in order to connect affine Gaudin models with integrable field theories we must now deviate slightly from the finite setting. Indeed, the Lax matrix L (z, x) of a classical integrable field theory is naturally thought of as the spatial component of a connection ∂ x + L (z, x) on the circle. For this reason, we say that an integrable field theory is a realisation of a local affine Gaudin model if its Lax matrix is of the form L (z, x) = ϕ(z) −1 Γ(z, x), (1.8) in which the g-valued fields such as J (i) (x) are expressed in terms of the observables of the theory, and if its twist function is given by ϕ(z). Moreover, the Hamiltonian of the field theory should be expressible as a linear combination of the Hamiltonians of the local affine Gaudin model, just as in the finite-dimensional setting. It is clear from the relation (1.8) and the common pole structure of the twist function and Gaudin Lax matrix that the zeroes of ϕ(z) coincide with the poles of L (z, x). In fact, it turns out that these zeroes, which we assume are all simple and call ζ i for i = 1, . . . , M , determine a particularly nice set of Hamiltonians for the local affine Gaudin model defined by the densities q i (x) := −res ζ i 1 2 ϕ(z) −1 κ(Γ(z, x), Γ(z, x))dz. (1.9a) These are to be compared with (1.4), used in the finite setting. In fact, it was shown in [6] that the zeroes of the twist function also play an important role in describing the higher local charges in integrable field theories with twist function, beyond the quadratic one given by (1.9a).
In this article we shall consider the class of integrable field theories which are realisations of affine Gaudin models associated with g, whose twist functions have simple zeroes and with Hamiltonian density given by a linear combination M i=1 i q i (x), (1.9b) for some parameters i , i = 1, . . . , M .
We are now ready to formulate the central problem addressed in the present paper. Consider two integrable field theories belonging to the above class. We let ϕ i (z) and Γ i (z, x) for i = 1, 2 denote their corresponding twist functions and Gaudin Lax matrices. We would like to construct a new integrable field theory depending on a parameter γ = 0 which couples these together in such a way that we recover the two uncoupled theories in the limit γ → 0. The above analogy with finite-dimensional Gaudin models serves as a guide for how to do this. Indeed, one can use the exact same construction as in (1.7) to define the Lax matrix of the affine Gaudin model corresponding to this new integrable field theory. And from this one can then also extract the twist function ϕ(z) and Gaudin Lax matrix Γ(z, x) of this coupled field theory. However, since we are now working with a Hamiltonian given by a density of the form (1.9b) associated with zeroes of the twist function, as opposed to the Hamiltonians naturally associated with the poles u i as in (1.4), the coupling procedure in the affine case is more subtle. The details of the latter are formulated in Theorem 2.1, which is the first main result of the paper.
As already emphasised earlier, the description of integrable field theories in terms of affine Gaudin models is intrinsically Hamiltonian. It is therefore not clear, a priori, how to determine whether such a field theory is relativistic from its Gaudin model description. In particular, in the coupling procedure outlined above, even if one starts with two relativistic integrable field theories it is not obvious whether the resulting coupled field theory is also relativistic. Remarkably, for the class of integrable field theories under consideration, namely with Hamiltonian densities of the form (1.9b), there exists simple conditions ensuring Lorentz invariance. Specifically, we show that an integrable field theory within this class is relativistic if i = ±1 for each i = 1, . . . , M . This is the second main result of the paper.
The plan of this article is the following. In Section 2 we begin by recalling the general formalism of local affine Gaudin models from [4]. In paragraph 2.3.3, we then give a general discussion of how two affine Gaudin models can be coupled together to form a third affine Gaudin model depending on an extra parameter γ in such a way that the decoupling limit corresponds to γ → 0. We go on to discuss the space-time symmetries of such models in Section 2.4 and in particular derive a simple condition on the form of the Hamiltonian for the theory to be Lorentz invariant. In Section 3 we apply this formalism to construct the action presented in [1] for N coupled principal chiral models with Wess-Zumino terms. As a warm-up, we begin by recalling in Subsection 3.2 how the principal chiral model with a Wess-Zumino term, corresponding to the case N = 1, is described in the language of affine Gaudin models, before moving on to the case of general N in Subsection 3.3 where we also discuss the symmetries of the resulting model. In Section 4 we give the affine Gaudin model descriptions of various integrable σ-models that may be used as building blocks for constructing more general coupled integrable σ-models following the procedure presented in Section 2. Finally, we collate some technical results in five appendices. Specifically, Appendix A concerns the pole structure of the generating function of the quadratic Hamiltonians of a local affine Gaudin model. Appendix B contains two elementary lemmas on the interpolation of rational functions. In Appendix C we provide details of the coupling procedure presented in Section 2, including a proof of Theorem 2.1. In Appendix D we apply the general construction to couple together N − 1 copies of the principal chiral model each with a Wess-Zumino term and one homogeneous Yang-Baxter deformation of the principal chiral model. This can be seen as the result of taking the model constructed in Section 3 with the Wess-Zumino term removed for one of the principal chiral fields, and then applying a homogeneous Yang-Baxter deformation to the latter. Finally, Appendix E contains a summary of Section 4 in the form of a table, listing various integrable σ-models and their affine Gaudin model descriptions that can be used as elementary building blocks in the general construction of Section 2.

Realisations of local Affine Gaudin models
This section is devoted to the study of Affine Gaudin Models (AGM) and their realisations. The first two subsections consist mainly in a review of the article [4] (see also [7] for a detailed review). In Subsection 2.1, we shall in particular discuss Takiff currents and realisations of Takiff algebras. In Subsection 2.2 we start by recalling the definitions of the Gaudin Lax matrix and the twist function and their relation to the usual Lax matrix which satisfies a Maillet non-ultralocal bracket. We also review the construction of quadratic Hamiltonians from which the Hamiltonian is extracted. The dynamics governed by this Hamiltonian takes the form of the zero curvature equation for a Lax pair, whose expression is given. We end this Subsection by explaining the presence of a global diagonal symmetry for all local AGM. In Subsection 2.3, we describe the parameters characterising a local AGM and its realisations. We show then how it is possible to couple two realisations of AGM and indicate the mechanism which enables to recover the decoupled limit. The Subsection 2.4 is devoted to the study of space-time symmetries. Since the construction of local AGM is at the Hamiltonian level, the analysis carried out in that Subsection allows to choose the parameters of a local AGM in order to obtain a field theory which is relativistic invariant.

Definitions and conventions
Lie algebras. Let g be a finite-dimensional simple complex Lie algebra. We denote by κ the opposite of its Killing form, which is a non-degenerate bilinear form on g. Let (I a ) a∈{1,··· ,n} be a basis of g and (I a ) a∈{1,··· ,n} be its dual basis with respect to κ. We then consider the split quadratic Casimir of g, defined as the element C 12 = I a ⊗ I a (2.1) of g ⊗ g, which does not depend on the choice of basis (I a ) a=1,··· ,n . The ad-invariance of the form κ then translates to the following identity on C 12 , where we use the standard tensorial index notations: Let us also note the following completeness relation: Let us fix a real form g 0 of g. It can be seen as the subalgebra of fixed point g 0 = {X ∈ g, τ (X) = X} of an antilinear involutive automorphism τ of g. We will suppose that the basis (I a ) a∈{1,··· ,n} of g chosen above is such that all I a 's belong to g 0 (there always exists such a basis): it then forms a basis of g 0 over R. The split quadratic Casimir (2.1) of the algebra is real, in the sense that it satisfies Finally, let us note that if we choose g 0 to be the compact form of g, the bilinear form κ reduces to a positive scalar product on g 0 .
Hamiltonian field theory. In this article, we will consider field theories on a one-dimensional space D, which will be either the real line R, parametrised by a coordinate x in ] − ∞, +∞[, or the circle S 1 , parametrised by a periodic coordinate x in [0, 2π[. We denote by δ xy = δ(x − y) the Dirac distribution on D and by δ xy = ∂ x δ(x − y) its derivative. For any field φ(x), we then have D dy φ(y)δ xy = φ(x) and D dy φ(y)δ xy = ∂ x φ(x).
We will consider Poisson algebras describing the local observables of Hamiltonian field theories on D. Such an algebra A contains fundamental dynamical fields φ 1 (x), · · · , φ k (x): A is then generated by these fields, in the sense that it contains algebraic combinations of these fields and their derivatives (evaluated at the same point x ∈ D) and integrals of such combinations over D 5 . More physically, this means that A contains the local fields and local charges (integrals) constructed from the fundamental fields φ i (x). The Poisson structure on A is then entirely determined by the Poisson bracket between the fundamental fields φ i , as one can then extend the bracket to the whole algebra A by Leibniz rule and linearity.
Suppose we are given two such Poisson algebras A and A, with fundamental fields φ 1 (x), · · · , φ k (x) and φ 1 (x), · · · , φ l (x). One can consider the algebra generated (in the sense above) by the fields φ i and φ i , equipped with a Poisson structure which extends the ones between the φ i 's and between the φ i 's by requiring Mathematically, this algebra is then the tensor product A ⊗ A (or more precisely an appropriate completion of this tensor product). In the rest of this article, we will understand tensor products between Poisson algebras in this sense.

Phase space: Takiff currents
Sites and levels. The phase space of a local AGM is characterised by the following data: • two finite abstract sets Σ r and Σ c , whose elements we call the real and complex sites of the model, for reasons to be explained below ; • for any α ∈ Σ r Σ c , a positive integer m α ∈ Z ≥1 , that we call multiplicity of the site α ; • for any α ∈ Σ r Σ c , a collection of m α numbers α p , p ∈ {0, · · · , m α − 1}, called the levels of the site α, that belong to R if α ∈ Σ r and to C if α ∈ Σ c .
For reasons to be explained in Subsubsection 2.1.3, we shall also suppose that for all sites α, the highest levels α mα−1 are non-zero. For compactness, we will gather these data in a unique object: that we will call a Takiff datum. For any complex site α ∈ Σ c , we introduce a labelᾱ, which we call the conjugate site to α. We then form the set Σ c of conjugate sites and the set Σ = Σ r Σ c Σ c of all sites. For all conjugate sitesᾱ ∈ Σ c , we define the corresponding levels ᾱ p = α p , for p ∈ {0, · · · , mᾱ − 1}, with multiplicity mᾱ = m α .
Takiff currents and observables. With each site α ∈ Σ, we associate m α currents J α [p] (x), p ∈ {0, · · · , m α − 1}. If α is a real site in Σ r , we suppose that these currents are real in the sense that they belong to g 0 (and are thus invariant under τ ). Moreover, we suppose that the currents associated with a complex site α in Σ c are related to the ones associated with the conjugate siteᾱ ∈ Σ c by τ J α [p] = Jᾱ [p] . The algebra of observables T of the local AGM is then defined as the algebra generated by the components of all currents J α [p] (in the sense defined in Paragraph 2.1.1). One then defines a Poisson structure on T by specifying the Poisson bracket between these currents: for all α, β ∈ Σ, p ∈ {0, · · · , m α − 1} and q ∈ {0, · · · , m β − 1}.
These are called Kac-Moody currents.

Generalised Segal-Sugawara integrals and momentum
Generalised Segal-Sugawara integrals. Let us fix a site α ∈ Σ. Following [4], we define complex numbers η α p for p ∈ {0, · · · , 2m α − 2} by the set of linear equations (see Appendix A, Lemma A.2 for more details about the numbers η α p ): As explained in [4], Lemma 4.1, these equations admit a unique solution (η α p ) p∈{0,··· ,2mα−2} , which can be expressed in terms of the levels ( α p ) p∈{0,··· ,mα−1} . It is this result which necessitates that the highest levels α mα−1 are non-zero for all sites α ∈ Σ. Moreover, one has (see Appendix A, Remark 1) For p ∈ {0, · · · , m α − 1}, we then define the quantity that we call generalised Segal-Sugawara integral (note that the sum on q and r can actually be restricted to the domain q + r ≥ m α − 1 + p as η α s = 0 for 0 ≤ s ≤ m α − 2). As a local charge constructed from the fundamental fields J α [p] (x) of the model, D α [p] is then an element of the algebra of observables T of the Gaudin model. Starting from the Takiff Poisson algebra (2.3), one then shows [4] that is complex and its conjugate is given by D α [p] = Dᾱ [p] , for all p ∈ {0, · · · , m α − 1}.
Cases of multiplicity 1 and 2.
For the examples of local AGM that we shall consider in this article, we will consider sites with multiplicities at most 2. Let us then describe more explicitly the construction of the Segal-Sugawara integrals D α [p] for multiplicities m α equal to 1 and 2. Let us start with the case of multiplicity m α = 1. In this case, one only has one level α 0 and one Takiff current J α [0] (which is then a Kac-Moody current). We are then searching for only one number η α 0 , which satisfies η α 0 α 0 = 1 (note that this equation admits a solution if the highest level α 0 is non-zero, as expected from the general discussion above). The corresponding Segal-Sugawara integral is then given by The density of this charge is a classical version of the Segal-Sugawara energy-momentum tensor associated with the Kac-Moody current J α [0] (x), which is used in the context of Conformal Field Theory (see for instance [8]). This observation justifies the name of generalised Segal-Sugawara integrals for the case with an arbitrary multiplicity.
Let us now investigate the case of multiplicity m α = 2 (which can be seen as the classical limit of the construction [9]), for which we have two levels α 0 and α 1 . One then has to search for three numbers η α p , for p = 0, 1, 2, satisfying One easily checks that this system of equations has a unique solution, given by Note that this requires that the highest level α 1 is non-zero and that η α 0 = 0, as expected from equation (2.5). One then easily constructs the two Segal-Sugawara integrals by applying Equation (2.6): Momentum of the local AGM. Let us consider the sum One easily checks from the Poisson bracket (2.7) that As all fields in the algebra T are combinations of the Takiff currents J α [p] (x), the Hamiltonian flow of P then generates the spatial derivative with respect to x on all these fields. Hence, P is the momentum of the local AGM. Note from the second expression in Equation (2.10) that P is real, as

Realisations of local Affine Gaudin Models
Realisations of Takiff algebra. Let us consider a Poisson algebra A, describing the observables of a field theory on D. We say that A is a realisation of T if there exists a (conjugacy-equivariant) Poisson map π : T −→ A.
Concretely, finding such a realisation of A amounts to constructing commuting Takiff currents from the dynamical fields in A with appropriate reality conditions. More precisely, one has to construct currents J α [p] (x) for all α ∈ Σ and p ∈ {0, · · · , m α − 1}, satisfying the same Poisson bracket as the abstract Takiff current J α [p] (x), i.e. such that where {·, ·} A denotes the Poisson bracket on A. In addition, we ask for the currents J α [p] to satisfy the same reality conditions as the abstract Takiff currents and thus require that The map π defining the realisation is then given by , ∀ α ∈ Σ and ∀ p ∈ {0, · · · , m α − 1}.
As the currents J α [p] generate the algebra T , these relations define entirely the map π, which is then a Poisson map as the J α [p] satisfy the same Poisson bracket as the J α [p] . The condition (2.13) then ensures that the map π is conjugacy-equivariant, i.e. that π(X) = π X , for any observable X ∈ T .
Suitable realisations. As we saw in Paragraph 2.1.3, the momentum of the local AGM is given by P . Applying the Poisson map π to the bracket (2.11), we get Thus, the Hamiltonian flow of π P coincides with the one of the momentum P A of A on the currents J α [p] (x). We will say that the realisation π is a suitable realisation if π P = P A . (2.14) Note that this condition is not automatically satisfied by all realisations, as the algebra A can contain more fields than the currents J α [p] (x), on which the Hamiltonian flow of π P may not generate the spatial derivative.
Physical models as realisations of local AGM. The notion of realisation presented above allows one to realise the observables T of a local AGM in terms of the observables of another algebra A. The reference [4] studies the local AGM on T and constructs for it an integrable dynamics, in particular a Hamiltonian whose equations of motion can be recast as a zero curvature equation, with a Lax matrix satisfying a Maillet bracket. All these constructions can then be transferred to the realisation A, hence allowing the construction of an integrable field theory on A.
In the Subsection 2.2, we will explain in details the construction of such an integrable model. For simplicity, we will consider directly the construction in A, i.e. in the realisation of T . Note however that all this construction also applies to the abstract local AGM on T .
Considering realisations of local AGM instead of the local AGM itself allows the description of a larger class of integrable models (although these can be seen as "images" of AGM). In particular, this allows to construct models which also possess a Lagrangian formulation. Indeed, the abstract Takiff algebra T cannot be understood as the algebra of Hamiltonian observables of a Lagrangian field theory, generated by coordinates fields φ i (x) and their conjugate momenta π i (x). However, it is possible in some cases to find a realisation A of T which comes from such a Lagrangian field theory. In this case, one can then perform the inverse Legendre transform and write an action for the considered model. This will be for example the case for the integrable σ-models that we will consider in Sections 3 and 4 of this article.
Finally, let us note that the notion of suitable realisation introduced in the previous paragraph will allow us to discuss space-time symmetries of realisations of local AGM in Subsection 2.4. In particular, as a suitable realisation provides a reinterpretation of the momentum P A of the theory in terms of the Takiff currents J α [p] (x), this will allow us to give a simple condition on the parameters defining the integrable model on A for this model to be Lorentz invariant.
For complex sites α ∈ Σ c , we also associate a position zᾱ = z α with the conjugate siteᾱ ∈ Σ c . By a slight abuse of language and as it is commonly done in the literature about Gaudin models, we will sometimes refer to the positions of the sites as the sites themselves.
Gaudin Lax matrix and twist function. One of the fundamental ingredients of a local AGM is its Gaudin Lax matrix. We define it here for the realisation A of the local AGM as It is a g-valued field in A, depending on an auxiliary parameter z ∈ C, called the spectral parameter.
Note the analytical structure of Γ(z, x) as a rational function of the spectral parameter: it possesses poles at the positions z α of the sites, whose orders are the corresponding multiplicities m α . Let us also introduce the following function of the spectral parameter: where ∞ is a real number that we will suppose to be non-zero 6 . This function ϕ(z) is the so-called twist function of the model. Except for the constant term ∞ , its analytical structure is similar to the one of Γ(z, x). One easily checks from the reality conditions on the currents and the levels that the Gaudin Lax matrix and the twist function obey the following (conjugacy-equivariance) reality conditions: τ Γ(z, x) = Γ(z, x) and ϕ(z) = ϕ(z). (2.17) As explained in [4], the Gaudin Lax matrix Γ(z, x) and the twist function ϕ(z) are both parts of the formal Lax matrix of the AGM (which is not to be confused with the usual Lax matrix of a field theory, as in Equation (2.19) below). An important characteristic of this formal Lax matrix (see [4], Proposition 4.12) is that it satisfies a linear Skyanin bracket in the affine untwisted Kac-Moody algebra associated with g. We shall not enter into this formalism in this article and will just use the fact that this bracket translates as the following bracket on the Gaudin Lax matrix: as can be also shown directly from the Poisson brackets (2.12) (for simplicity, we dropped the subscript A for the Poisson bracket on A as we are now working exclusively in this fixed realisation).
Lax matrix, Maillet bracket and R-matrix. We define the Lax matrix of the model as As Γ(z, x), it is a g-valued field in A, depending rationally on the spectral parameter and satisfying the reality condition: τ L(z, x) = L(z, x).
The Poisson bracket (2.18) on Γ(z, x) translates into the following bracket on L(z, x): The bracket (2.20) is a Maillet non-ultralocal bracket, as introduced in [10,11]. The fact that it satisfies the Jacobi identity is ensured by the fact that the R-matrix is a solution of the classical Yang-Baxter equation We end this paragraph by the following remark, which will be useful for Section 3. From Equations (2.15), (2.16) and (2.19), one gets for all α ∈ Σ.

Quadratic Hamiltonians and zero curvature equation
Quadratic Hamiltonians. Let us define the quantity (2.24) are the generalised Segal-Sugawara integrals in the realisation A. The quantity P(z) has the same analytical structure than the Gaudin Lax matrix Γ(z, x) and the twist function ϕ(z): as explained in [4], this comes from the fact that together, they form the formal Lax matrix of the model. It satisfies the reality condition P(z) = P(z).
From the bracket (2.7) (and applying the Poisson map π), one finds that In [4], Propositions 4.14 and 4.15, it is explained how to obtain quadratic charges in involution for local AGM as part of a spectral parameter dependent Hamiltonian constructed from the formal Lax matrix. Here, we will express it, for the realisation A, in terms of the Gaudin Lax matrix Γ(z, x), the twist function ϕ(z) and the quantity P(z) as It is a local charge in the algebra of observables A, depending rationally on the spectral parameter, whose density is a quadratic polynomial in the Takiff  for all z, w ∈ C. This construction then gives a set of quadratic Hamiltonians in involution as the quantities extracted linearly from H(z): evaluations at special values of z, residues (or more generally coefficients of poles in z), combinations of those, ...
A particular way of extracting a complete set of such Hamiltonians is to determine the partial fraction decomposition of H(z) seen as a rational function of z. It is clear from the expressions (2.15), (2.16) and (2.23) of Γ(z), ϕ(z) and P(z) that H(z) has poles only at the positions z α of the sites α ∈ Σ and of order at most 2m α , with m α the multiplicity of the site α. An important result for Subsection 2.2.3 is that these poles are actually of order at most m α . This is a consequence of the definition of the generalised Segal-Sugawara integrals in Paragraph 2.1.3. As the demonstration of this result is rather technical, we detail it in the Appendix A. As a consequence, we can write the quadratic Hamiltonian H(z) as (2.28) We then have H α [p] , H β [q] = 0, for all α, β ∈ Σ, p ∈ {0, · · · , m α − 1} and q ∈ {0, · · · , m β − 1}. Moreover, any quantity linearly extracted from H(z) can be seen as a linear combination of the H α [p] 's. In particular, we shall choose the Hamiltonian of the model as such a linear combination: for some γ α p ∈ C. The dynamics of the model is then given by the Hamiltonian flow   We shall not enter into more details about the expression of the M α [p] (z, x)'s, as in Subsection 2.2.3 we will give another parametrisation of the Hamiltonian yielding a simpler expression of M(z, x).
Integrability. The fact that the equations of motion of the model can be recast as a zero curvature equation (2.32) implies the existence of an infinite number of conserved quantities for this model, extracted from the monodromy matrix of the Lax matrix L(z, x). Moreover, the fact that the Lax matrix satisfies the Maillet bracket (2.20) ensures that these conserved charges are pairwise in involution. Thus, realisations of local AGM are by construction integrable two dimensional field theories.
As the R-matrix governing the Maillet bracket is of the form (2.21), they in fact belong to a more restrictive class of integrable field theories, called models with twist function [12,13,14,15] (see also [7]). This additional structure allows a deeper understanding of these models and their properties. In particular, one can apply the method developed in [6] to construct an infinite integrable hierarchy of local charges. We will come back to this below.

Zeros of the twist function and integrable hierarchy
Zeros of the twist function. Let us consider the twist function (2.16) of the model. One can write it as the quotient of two co-prime polynomials of degree as we assumed that ∞ is non-zero. Thus, the twist function admits M zeros ζ i , i ∈ {1, · · · , M }, in the complex plane. One can then write For the rest of this article, we will suppose that the zeros ζ i of ϕ(z) are simple.
Extracting Hamiltonians from the zeros. Recall the spectral parameter dependent Hamiltonian H(z) introduced in (2.26). As Γ(z) and P(z) are regular at ζ i (because the zeros ζ i of the twist function are all different than its poles z α ), the evaluation of H(z) at ζ i is given by As H(z) possesses the partial fraction decomposition (2.28), it is clear that the H(ζ i )'s, i ∈ {1, · · · , M }, are linear combinations of the H α [p] 's, α ∈ Σ and p ∈ {0, · · · , m α − 1}. As we supposed that the twist function has simple zeros, so that the ζ i 's are pairwise distinct, the converse is also true according to the auxiliary Lemma B.1. Thus, the sets of Hamiltonians H(ζ i ) i∈{1,··· ,M } and H α Instead of the evaluations H(ζ i ), we shall use slightly different quantities in the rest of this article, constructed as follows. Let us introduce the local charge We then define for all i ∈ {1, · · · , M }. Recall that we suppose the zeros ζ i of ϕ(z) to be simple. Then, as Γ(z) is regular at ζ i , we have where ϕ (z) denotes the derivative of the twist function with respect to the spectral parameter (note that ϕ (ζ i ) is non-zero as ζ i is a simple zero of ϕ(z)). The charges Q i are then pairwise in involution: The set of charges (Q i ) i∈{1,··· ,M } forms a complete set of quantities extracted linearly from H(z). In particular, the Hamiltonian H defined in Equation (2.29) can be written as a linear combination for some numbers i . For the rest of this article, we shall use this parametrisation of the Hamiltonian of the model, rather than the one (2.29).
Lax pair and zeros. Recall the Lax matrix L(z, x) of the model, defined by Equation (2.19). As ϕ(z) and Γ(z, x) have poles at the same positions z α , α ∈ Σ, and with the same orders m α (as we required that highest levels α mα−1 are non-zero), the Lax matrix has no singularities at the z α 's. Thus, it has only singularities at the zeros of the twist function. Moreover, as the zeros of the twist function are simple, these singularities are simple poles. The Lax matrix can then be written Let us now consider the temporal component of the Lax pair, M(z, x). It can be expressed as (2.33) if one uses the parametrisation (2.29) of the Hamiltonian H. As explained in this subsection, we shall now use the parametrisation (2.38) of H instead, using the zeros of the twist function. From Equation (2.37) and the spectral parameter dependent zero curvature equation (2.30), we see that the matrix M(z, x) can be written in this parametrisation as From the expression (2.31) of M(w ; z, x), we then deduce that We note the similarity of this expression with the one (2.39) of the Lax matrix. This is one of the motivations for the use of the parametrisation (2.38) of the Hamiltonian, as it allows to treat the spatial and temporal components of the Lax pair in similar ways. This and related ideas will be useful to treat the space-time symmetries and Lorentz invariance of AGM in Subsection 2.4.
Positivity of the Hamiltonian. For this paragraph only we shall suppose that the real form g 0 of the Lie algebra g is the compact form. The bilinear form κ restricted to g 0 is then positive. Let us also suppose that the zeros ζ i of the twist function are real. We then deduce from the reality conditions (2.17) that Γ(ζ i , x) belongs to the real form g 0 and ϕ (ζ i ) is real. In particular, the evaluation H(ζ i ), given by (2.34), is positive for all configurations of the fields and the charge Q i , given by (2.37), has the sign of −ϕ (ζ i ). Thus, if we choose i to be of this same sign, the Hamiltonian H is positive (and in particular bounded below).
Local integrable hierarchy. Let us end this subsection with a brief discussion about integrable hierarchies in realisations of local AGM. It was shown in [6] that models with twist function possess infinite integrable hierarchies of local charges associated with the zeros of the twist function 7 . These charges are given by where the Φ d 's are well-chosen invariant polynomials [16] on g of degree d + 1. The integers d are restricted to specific values, namely the exponents of the untwisted affine algebra of g. In particular, all affine algebras possess d = 1 as an exponent, associated with quadratic charges Q 1 i . The corresponding polynomial Φ 1 is then the bilinear form κ. Thus, these quadratic charges Q 1 i coincide with the Hamiltonians Q i introduced above.
The charges Q d i , for i ∈ {1, · · · , M } and d running over the affine exponents of g, are pairwise in involution. As the Hamiltonian H of the model is a linear combination of the quadratic charges Q 1 i = Q i , these charges are also conserved. We then get M infinite towers of local conserved charges in involution. Moreover, it is shown in [6] that the Hamiltonian flows of these charges on the Lax matrix L(z, x) take the form of compatible zero curvature equations, generalising the result of Paragraph 2.2.2 for quadratic Hamiltonians. This forms what is called an integrable hierarchy.

Global diagonal symmetry of local AGM
Conserved charge. Let us consider the field By the reality condition (2.13), it is valued in the real form g 0 . From Equation (2.30), one finds that K ∞ is the time component of a conserved current, as and as the Hamiltonian of the model is extracted linearly from H(z). Thus, the g 0 -valued local charge is a conserved charge of the model. Moreover, from Equation (2.18), one shows that K ∞ satisfies the Kirillov-Kostant bracket of g 0 : Associated diagonal symmetry. By the Noether theorem, this conserved charge is associated with an infinitesimal symmetry of the model with underlying symmetry algebra g 0 . The corresponding variation of an observable O ∈ A with infinitesimal parameter ∈ g 0 is given by This transformation leaves invariant the Hamiltonian H of the model and more generally the spectral parameter dependent Hamiltonian H(z), as seen by integrating equation (2.42) over x: 7 In [6], hierarchies are associated more precisely with so-called regular zeros ζi of the twist function, for which the Gaudin Lax matrix Γ(z, x) is regular at z = ζi. In the context of local AGM and their realisations, all zeros of the twist function are automatically regular, as they are different from the positions zα of the sites α ∈ Σ. Note also that the result of [6] was proved only for algebras g of classical types A, B, C and D (although we do believe that a similar result should exist for exceptional types).
From Equation (2.18), one shows that the variation of the Gaudin Lax matrix is given by As a consequence, one has, using Equation (2.15), We call δ ∞ the diagonal symmetry of the model, as it acts similarly on each Takiff current J α [p] . This infinitesimal symmetry lifts to a global action of the adjoint group G 0 of g 0 : By invariance of the form κ, one then sees that the generalised Segal-Sugawara integrals D α [p] are also left unchanged under this transformation. Similarly, by invariance of the polynomials Φ d , one gets that the whole integrable hierarchy of local charges Q d i defined in Equation (2.41) is invariant under this transformation, generalising the invariance of the quadratic charges Q 1 i .

Parameters of the model
Let us summarise what are the parameters characterising the model introduced in the two previous subsections: • the Takiff datum (which contains the data of the sites, their multiplicities and their levels), defining, together with the Takiff realisation π, the algebra of observables of the model ; • the positions z of the sites ; • the constant term ∞ in the twist function ; • the coefficients = ( i ) i=1,··· ,M defining the Hamiltonian through Equation (2.38).
Note that the datum of , z and ∞ is exactly equivalent to the choice of a twist function ϕ(z) as in (2.16). To stress this dependence, we will sometimes denote the Hamiltonian (2.38) as H = H ϕ,π when we have to discuss several (realisations of) affine Gaudin models at the same time. Similarly, we will sometimes use the notation M ϕ,π to identify in a compact way the model we are considering. Note that there is a constraint on the parameters defining M ϕ,π . Indeed, to describe the Hamiltonian of the theory from the quadratic charges extracted from the zeros of the twist function, we supposed that these zeros are simple.

Change of spectral parameter
Change of coordinates in the spectral plane. A crucial element in the construction of local AGM and their realisations is the spectral parameter z, which plays the role of an arbitrary auxiliary parameter. It is then natural to wonder if there are other possible choices of spectral parameter which preserve the structure of local AGM.
Let us first determine what would be such a possible change of spectral parameter. It would be mathematically represented as an invertible map from the complex plane C to itself. As we are considering rational functions of the spectral parameter, this map should send rational functions to rational functions and thus be rational itself. Yet, an invertible rational map from C to itself is a linear polynomial, i.e. an affine transformation 8 . Moreover, this map should send the positions of the sites z to another set of admissible positions: in particular, it should send real points to real points and pair of conjugated points to pair of conjugated points. We shall then consider the maps with a ∈ R * and b ∈ R.
Transformation of the twist function. Let us then fix such a map and consider the change of spectral parameter z → z = τ a,b (z) = az + b. In the construction of local AGM presented in [4], the twist function appears naturally as a differential 1-form ϕ(z) dz. Thus, the twist function should transform as Considering the expression of the affine map τ a,b , one simply gets here that Let us consider the partial fraction decomposition (2.16) of the twist function ϕ(z). One then finds from equation (2.43) a similar decomposition for the new twist function: Transformation of the realisation. The model initially considered is defined from a realisation π : T → A of the Takiff algebra with levels . This realisation was given by Takiff currents J α [p] (x) in A (see Paragraph 2.1.4). It is easy to check that the currents (2.46) satisfy the Poisson brackets of Takiff currents with levels . This then defines a realisation π : T −→ A of the algebra with Takiff datum . It is easy to check that the the corresponding Gaudin Lax matrix is given by Thus, the Gaudin Lax matrix also obeys the transformation law of a 1-form: which is coherent, as it also appears naturally as a 1-form in the reference [4].
Transformation of the quadratic Hamiltonians. Let us consider the images ζ i = τ a,b (ζ i ), i ∈ {1, · · · , M } of the zeros ζ i of the twist function ϕ(z). By construction, they are the zeros of the transformed twist function ϕ( z). Let us also consider the quantity Q(z) defined in (2.35). Considering the transformations above of the twist function and the Gaudin Lax matrix, one sees that this quantity transforms according to and thus as a 1-form, i.e. such that be the equivalent for the transformed model of the quadratic charges Q i defined in (2.36). As residues of 1-form are invariant under change of coordinates, the transformed charges Q i actually coincide with the initial charges Q i . In particular, we then have the equality of the Hamiltonians of the models: without redefining the parameters . Note that other sets of quadratic charges, such as the evaluations or the residues of H(z), would not be invariant under the change of spectral parameter. This also motivates the choice of the Q i 's as a "good" basis of quadratic charges and of as "good" parameters of the model.
after the change of spectral parameter.
Summary. As a conclusion, we have proved that the models M ϕ,π and M ϕ, π are identical, as they possess the same Hamiltonian. This then exhibits a redundancy in the set of parameters ϕ, π, used to describe the model. Let us discuss this redundancy in more details. Dilations τ a,0 modify the levels of the model through the relation (2.45). Thus, they modify the data defining the algebra of observables of the local AGM. However, this modification is counterbalanced by the one of the realisation π through Equation (2.46), so that the two models still share the same algebra of observables A. Moreover, as we see in Equation (2.44), dilations also modify the level at infinity ∞ , i.e. the constant term in the twist function, by a factor 1/a. One can then use the redundancy coming from dilations to fix the value of ∞ .
Translations τ 0,b act differently on the parameters of the model. Indeed, they only modify the positions z α of the sites. Thus, they can be used to fix the position of one site to a particular point in the complex plane, without changing the Takiff datum of the model.

Coupling and decoupling realisations of local AGM
Combining two Takiff realisations. Let us consider two Takiff data for k ∈ {1, 2}. We denote the levels of both Takiff data using the same notation α p , the distinction between the two being contained in the fact that the labels of the sites α are taken in different sets Σ c and consider the Takiff datum which contains the sites and levels of both 1 and 2 . The Takiff algebra T 1⊗2 is generated by Takiff currents J α [p] (x) for α in Σ (1) and Σ (2) . The Poisson bracket between the Takiff currents for α ∈ Σ (k) is the same as in the algebra T k of k alone. Moreover, the Poisson bracket between Takiff currents of sites in Σ (1) and Σ (2) is zero. Thus, the algebra T 1⊗2 can be seen as the tensor product as defined in Paragraph 2.1.1. Let us suppose that we have realisations π k : T k → A k for both k equal to 1 and 2. Then, the map defines a realisation. Moreover, this realisation is suitable (see Paragraph 2.1.4) if and only if both realisations π 1 and π 2 are.
Coupling realisations of AGM. Let us consider two models M ϕ 1 ,π 1 1 and M ϕ 2 ,π 2 2 , with corresponding Takiff data 1 and 2 (and constructed from the same underlying Lie algebra g). They are integrable models, with observables A 1 and A 2 respectively. In this subsection, we aim to construct an integrable model coupling these two models. The algebra of observables of such a model would be the algebra generated by the observables in both A 1 and A 2 , i.e. the tensor product A 1 ⊗ A 2 . As seen above, this tensor product is a realisation of the Takiff algebra T 1⊗2 , through the map π 1⊗2 . We shall then construct the coupled model as a realisation of a local AGM with Takiff datum 1⊗2 and sites Σ (1⊗2) = Σ (1) Σ (2) .
According to Subsection 2.3.1, the next step in the construction of the coupled model is to choose the positions of the sites and the level at infinity ∞ , which together with the Takiff datum 1⊗2 will determine the twist function of the model. At first sight, it could seem difficult to choose in a canonical way the constant term ∞ , as one could for example take the one of the twist function ϕ 1 or the one of the twist function ϕ 2 . However, this difficulty can be overcome using dilations of the spectral parameter. Indeed, as explained in Subsection 2.3.2, such a dilation for the model M ϕ 2 ,π 2 2 would modify the constant term in ϕ 2 without changing the model (and in particular its Hamiltonian). Thus, one can perform an appropriate dilation so that this constant term agrees with the one of ϕ 1 . We shall then choose ∞ to be this common value. Note that this dilation changes the other levels belongs to Σ (1) or Σ (2) . Remember however from Subsection 2.3.2 that one can modify these positions without changing the models M ϕ k ,π k k , by applying translations to their spectral parameters. Thus, one can define the positions of the sites of the coupled model to be where a 1 and a 2 are two real numbers. Note that one can also act by a translation of the spectral parameter on the coupled model too, thus shifting all positions w α by a common number. Hence, one can choose a 1 = 0. For future convenience, we shall parametrise a 2 as 1/γ, with γ a real number that we will call the coupling parameter. In conclusion, we have We will denote by ϕ 1⊗2,γ the twist function of the coupled model, thus defined by the Takiff datum 1⊗2 , the positions w α and the constant term ∞ . The last parameters to fix to define entirely the coupled model are the coefficients defining the Hamiltonian of the model. We will discuss how to choose these coefficients canonically in the next paragraph.
Coupled Hamiltonian. Let us consider the parameters 1 and 2 entering the definition of the models M ϕ k ,π k k . They are associated respectively with the zeros of the twist functions ϕ 1 and ϕ 2 . We define There are then M 1 zeros of ϕ 1 , that we shall denote by ζ We then label the parameters in 1 and 2 in a similar way: Let us note that to consider the models M ϕ k ,π k k , we have to suppose that the zeros ζ Let us note that there is a model on A 1 ⊗ A 2 which is also automatically integrable: the model with Hamiltonian H ϕ 1 ,π 1 1 + H ϕ 2 ,π 2 2 , corresponding to the two models M ϕ k ,π k k combined without any coupling. As we want to construct a model which couples M ϕ 1 ,π 1 1 and M ϕ 2 ,π 2 2 , we are looking for a choice of parameters such that the decoupled model is obtained in some limit of these parameters. The important result of this Subsection is then the following.
Moreover, for γ small enough, these zeros are simple. Besides, if we choose the parameters i to be equal to the corresponding Proof. As the proof of this theorem is rather long and technical, we give it in Appendix C. . It first ensures that, at least for small enough values of γ, the twist function ϕ 1⊗2,γ is simple and hence that the Hamiltonian of the corresponding realisation of local AGM is determined by parameters . Moreover, it provides a canonical choice of these parameters such that the limit γ → 0 of the model M ϕ 1⊗2,γ , π 1⊗2 corresponds to the decoupled combination of M ϕ 1 ,π 1 1 and M ϕ 2 ,π 2 2 . For γ = 0 but small enough, we thus constructed an integrable model on A 1 ⊗ A 2 , which couples non trivially M ϕ 1 ,π 1 1 and M ϕ 2 ,π 2 2 . This justifies a posteriori the name of coupling parameter for γ.
Coupled and decoupled Lax pair. Let us now discuss the Lax pair of the coupled model and its behaviour in the decoupling limit. By construction, the model M ϕ 1⊗2,γ , π 1⊗2 is integrable and possesses a Lax pair L γ (z), M γ (z) , satisfying a zero curvature equation in the Lie algebra g. Moreover, the decoupled model with Hamiltonian H ϕ 1 ,π 1 is also integrable and admits a natural Lax pair , for k equal to 1 or 2, is the Lax pair of the model M ϕ k ,π k k . It would seem at first that the Lax pairs of the coupled and decoupled model are unrelated, as they do not belong to the same Lie algebra. In fact, as shown in Appendix C, one has Similarly to what is discussed above (see Appendix C for more details), one has As a conclusion, we have the following.
form a Lax pair of the coupled model M ϕ 1⊗2,γ , π 1⊗2 . Moreover, in the decoupling limit, one has: (z) in Proposition 2.2 is necessary to get the decoupling limit (2.51), which cannot be obtained with only one of the g-valued factor.

Space-time symmetries
In this Subsection, we consider a realisation M ϕ,π of local AGM. From now on, we will suppose that the Takiff realisation π is suitable, as defined in Paragraph 2.1.4. This will be motivated by Proposition 2.3 which states that for such a suitable realisation, the momentum of the model has a simple expression in terms of the quadratic charges Q i . This will allow us to discuss various space-time symmetries of the model. In particular, we will give a simple characterisation of the relativistic invariance of the model in terms of the parameters of the model.

Space-time translations, energy-momentum tensor and scale invariance
Hamiltonian and momentum. Recall that the Hamiltonian of the model is defined as (2.38) in terms of the quadratic charges Q i . It does not depend explicitly on the time parameter t and is thus conserved. This is a sign of the invariance of the model under time translation. Moreover, the quadratic charge Q i is expressed as an integral of a local density, which depends on the space coordinate x only through the dynamical field Γ(ζ i , x) (see Equation (2.37)). Thus, Q i is invariant under space translation and so is the Hamiltonian H of the model. As a consequence, the momentum P A of the model, which generates the spatial derivatives of dynamical fields, is conserved. An important result of this Subsection is the following proposition, which gives a simple expression of the momentum P A in terms of the charges Q i . Proposition 2.3. We have (for a suitable realisation) Proof. Combining the definition (2.14) of a suitable realisation with Equations (2.10), (2.24) and (2.23), one gets From Equations (2.26) and (2.35), one gets Recall that H(z) has poles at the z α 's, of same orders m α as the twist function ϕ(z) (see Equation (2.28) and Appendix A). Thus the quantity H(z)/ϕ(z) has no residues at z = z α , for all α ∈ Σ. In particular, we get Recall the expression (2.15) and (2.16) of Γ(z, x) and ϕ(z). For small u, we have Thus, the 1-form is regular at u = 0, or in other words the 1-form Q(z) dz is regular at z = ∞. As a consequence, Q(z) dz has poles only at the positions z α of the sites α ∈ Σ and at the zeros ζ i of the twist function.
Recall that the sum of all residues of a 1-form on P 1 vanishes. Thus, we have Combining this with Equation (2.53), we get Energy-momentum tensor. For i ∈ {1, · · · , M }, we denote by q i (x) the density of the charge Q i . In this subsection, we will use greek labels µ, ν, ρ equal to 0 and 1 as space-time labels. In particular, we have Proposition 2.4. The components of the energy-momentum tensor of the model are Proof. The components T 0 0µ of the energy-momentum tensor are defined as the densities of the Hamiltonian and the momentum of the model: From Equations (2.38) and (2.52), we easily read As H and P A are conserved charges, their densities should satisfy local conservation equations, which serve as a definition of the components T 1 1ν : To obtain the expressions of these components, we will need the Poisson brackets of the densities q i . From Equation (2.7), one gets Thus, one has , this allows to read T 1 11 (x) as A similar computation yields Classical scale invariance. It is clear from Proposition 2.4 that the trace T µ µµ = T 0 00 + T 1 11 of the energy-momentum tensor vanishes. It is well known that this implies the classical scale invariance of the model. The conserved current associated with this invariance is The associated conserved charge 9 is

Relativistic invariance and light-cone formulation
Lorentz invariance. Let us consider the Minkoswki metric η = diag(+1, −1). It allows us to change contravariant indices into covariant ones and in particular to define It is clear from Proposition (2.4) that Thus, one sees that the model is relativistic if the parameters satisfy This remarkably simple condition is the main reason for the use of the parameters to define the Hamiltonian of the model. In the rest of this article, we shall assume that this condition is satisfied. We will need the following sets Light-cone coordinates and Lax pair. Let us introduce the light-cone coordinates We denote by ∂ ± = ∂ t ± ∂ x the corresponding derivatives. These derivatives can be seen as the Hamiltonian flows From the expressions (2.38) and (2.52) of H and P A , one gets (2.59) 9 More precisely, this charge is conserved only if we are considering a field theory on the real line, i.e. if D = R.
It is interesting to see that the charges Q i separate into two sets appearing respectively in P + for i ∈ I + and in P − for i ∈ I − . As all Q i 's Poisson commute, the charges P + and P − are in involution. One has in fact the stronger result that even the densities of P + and P − are in involution, as a consequence of Equation (2.56). The zero curvature equation (2.32) can also be rewritten in light-cone formulation as where we introduced the light-cone lax pair From the equations (2.39) and (2.40), we get Similarly to the expression (2.59), one sees that the expressions of L ± (z) separate the currents Γ(ζ i ) depending on whether i belongs to I + or I − . Moreover, they have a nice consequence on the analytic structure of the Lax pair in terms of the spectral parameter: although L(z) and M(z) possess poles at all ζ i 's for i ∈ {1, · · · , M }, the light-cone components L + (z) and L − (z) only have poles at the ζ i 's for i in I + and I − respectively. This will be an important feature for Section 3.
Spins of the currents. In this paragraph, we consider a model with space coordinate x on the real line R. Let F(x, t) be a local field (for this paragraph, we will write explicitly the dependence of the fields in t although we are still considering the model in the Hamiltonian formulation). Let us also consider a Lorentz boost with infinitesimal rapidity ϑ, which then changes the coordinates by δ B x = ϑ t and δ B t = ϑ x. Under this transformation, the field F(x, t) changes because of its dependence on the coordinates x and t but also because of its internal degree of spin. More precisely, we say that F(x, t) has spin s if its variation under the Lorentz boost is given by One then has the following proposition.
Proof. The transformation of Γ(ζ i , x, t) under the boost is given by where B is defined in (2.58). By Proposition 2.4, we get Thus, we have The first term in this equation is easy to express as {P A , ·} = ∂ x . The second term can be computed starting from the Poisson bracket (2.7) and the expression (2.54) of q i . For the sake of brievity, we will not detail this here. After a few manipulations, one finally obtains which ends the demonstration, as ϑ t = δ B x and ϑ x = δ B t.
Proposition 2.5 shows that the observables which are natural to consider to discuss space-time properties of the model are the currents Γ(ζ i ), as they possess a definite spin (which is equal to ±1 for i ∈ I ± ). This emphasises the importance of the zeros ζ i of the twist function, which was already observed in the previous paragraphs. Moreover, Proposition 2.5 proves that the local charges Q d i of the integrable hierarchy (2.41) have spin i (d + 1). Thus, the hierarchy attached to a zero ζ i with i ∈ I + (resp. i ∈ I − ) is composed of charges of positively increasing (resp. negatively decreasing) spins.
Coupling relativistic realisations of AGM. We end this section by a short discussion combining the ideas of the present subsection and the ones of Subsection 2.3.3 about coupled realisations of AGM. Let us then consider two suitable realisations of AGM M ϕ 1 ,π 1 1 and M ϕ 2 ,π 2 2 , following the notations of Subsection 2.3.3. We suppose that the parameters also yields a relativistic model.
Let us also discuss positivity properties of the Hamiltonians of these models. For that we restrict to the case of a compact real form g 0 and suppose that the zeros ζ to be equal to +1 or −1 for Lorentz invariance, we thus have As explained in Subsection 2.3.3, for γ small enough, the zeros ζ i (γ) of the coupled twist function ϕ 1⊗2,γ are in one-to-one correspondence with the ζ i 's to be real ensures, at least for γ small enough, that the ζ i (γ)'s are also real and that ϕ 1⊗2,γ ζ i (γ) and ϕ k ζ 3 New integrable σ-models 3.1 Fields on T * G 0 and Takiff realisation of multiplicity 2 Coordinates and conjugate momenta. Let us consider a real Lie group G 0 with Lie algebra g 0 , described by local coordinates 10 φ i for i ∈ {1, · · · , n = dim g 0 }. Let us also consider the cotangent bundle T * G 0 of G 0 . In addition of the coordinates φ i , it is locally described by momenta π i conjugate to these coordinates. Let us now consider a field on D valued in T * G 0 , given locally by scalar fields φ i (x) and π i (x). As T * G 0 is a cotangent bundle, it is equipped naturally with a symplectic structure. This transfers to the canonical Poisson bracket on the fields φ i (x) and π i (x): for all i, j ∈ {1, · · · , n}. We will denote by A G 0 the Poisson algebra generated by the fields φ i (x) and π i (x) (in the sense of Paragraph 2.1.1).
Group-valued and algebra-valued fields. We now want a coordinate-free description of this phase space. Acting by translation on the base G 0 of T * G 0 , one can always send each cotangent space T * p G 0 at p ∈ G 0 to the cotangent space T * Id G 0 at the identity Id ∈ G 0 . Yet, T * Id G 0 is nothing but the dual g * 0 of the Lie algebra g 0 . Moreover, as g 0 is a semi-simple Lie algebra, it admits the non-degenerate bilinear form κ, yielding a canonical isomorphism g * 0 g 0 . Hence, one has a canonical isomorphism between T * G 0 and G 0 × g 0 . Thus, one should be able to described the T * G 0 -valued field as a pair of field g(x) in G 0 and X(x) in g 0 . Let us describe more concretely these fields.
As the φ i 's are local coordinates on G 0 , the fields φ i (x) naturally recombine to form the G 0 -valued field g(x). Let us denote by ∂ i the derivative with respect to φ i . Then g −1 ∂ i g is a g 0 -valued object. Let us write it in a basis (I a ) a∈{1,··· ,n} of g 0 as The matrix L a i i,a=1,··· ,n is then invertible. We shall write L i a i,a=1,··· ,n its inverse. It verifies We then define the g 0 -valued field where (I a ) a=1,··· ,n is the dual basis of (I a ) a=1,··· ,n with respect to the bilinear form κ ab . One can check that this field X is independent of the choice of coordinates (φ i ) i∈{1,··· ,n} and basis (I a ) a∈{1,··· ,n} .
Coordinate-free Poisson structure. The canonical Poisson bracket (3.1) on the fields φ i (x) and π i (x) can be rewritten as a coordinate-free Poisson bracket on the fields g(x) and X(x). More precisely, we have Let us define the g 0 -valued current From the Poisson brackets (3.3), one finds that it satisfies the following brackets: Momentum. Let us define This is the momentum of the algebra of observables A G 0 . Indeed, from the Poisson brackets (3.3) and (3.4), one finds that its Hamiltonian flow generates the spatial derivative with respect to x on the fields g(x) and X(x): and

Wess-Zumino term
Canonical 3-form on G 0 . Let us recall the coordinates φ i on the group G 0 and the corresponding derivatives ∂ i . Let us introduce the following tensor: The ad-invariance of the bilinear form κ implies that the tensor ω is totally skew-symmetric. It defines a 3-form on G 0 by Ω = ω ijk dφ i ∧ dφ j ∧ dφ k .
One checks that this form is closed in the de-Rham cohomology of G 0 : Thus, it is locally exact and can be locally written as is a 2-form on an open subset of G 0 . In terms of coordinates, this translates to the following relation: The current W (x). Let us introduce the following g 0 -valued current: From the canonical Poisson bracket (3.1), one shows that this current satisfies and Note also that the skew-symmetry of λ ij implies that the currents j and W are orthogonal with respect to the Killing form κ: κ(j, W ) = 0. (3.10) Wess-Zumino term. For this paragraph, we will consider the fields to depend on a time coordinate t ∈ R in addition to the spatial coordinate x (as we are working in the Hamiltonian formulation, this time dependence is implicitly defined by the choice of an Hamiltonian). We will denote by W = R × D the corresponding space-time manifold parametrised by t ∈ R and x ∈ D. In particular, we can consider the G 0 -valued field g(x, t) on W. Let us extend it to a 3-dimensional manifold B whose boundary is ∂B = W, parametrised by the coordinates x and t along with a third coordinate ξ. We then get a field g(x, t, ξ). The so-called Wess-Zumino term is defined as [17,18,19] According to the definition (3.6) of the tensor ω ijk , this can be written in terms of the coordinates fields φ i as As the form Ω is closed, this 3-dimensional term on B can be written locally as a 2-dimensional term on the space-time W = ∂B: independent of the "bulk" coordinate ξ. From the definition (3.7) of the current W , one finally finds that Takiff currents. Let us fix two real numbers and (we shall suppose that is non-zero). We introduce the currents and J [1] (x) = j(x).
Moreover, these currents are valued in the real form g 0 and are thus real Takiff currents.
Suitable Takiff realisation. Let us consider a Takiff datum with a unique site, multiplicity two and levels as above. The currents introduced in the previous paragraph then define a Takiff realisation: from the corresponding Takiff algebra T to the algebra A G 0 of fields on T * G 0 . We shall call this the PCM+WZ realisation, for reasons to be made obvious in Subsection 3.2. In the particular case where the level is equal to zero, we shall speak of the PCM realisation. From the definitions (3.12) of the currents J [0] (x) and J [1] (x) and the orthogonality relation (3.10), one finds that We recognize on the right-hand sign of this equation the generalised Segal-Sugawara D [0] in the realisation A G 0 (comparing to the Equation (2.9) giving the expression of the generalised Segal-Sugawara integral of a Takiff current of multiplicity 2). Thus the PCM+WZ realisation π G 0 is suitable.

Warm-up: the Principal Chiral Model with Wess-Zumino term as a realisation of a one-site local AGM
In the previous subsection, we have exhibited a Takiff realisation in the algebra A G 0 of canonical fields on T * G 0 , the PCM+WZ realisation. Following the section 2 of this article, one can then construct an integrable field theory with observables A G 0 from the general construction of local AGM. In this subsection, we show that this integrable model can be identified as the Principal Chiral Model with Wess-Zumino term on the group G 0 , as initially observed in the article [4]. This will serve as a warmup for the rest of this section as it is the first and simplest example of the construction of integrable σ-models from AGM.

The model in Hamiltonian formulation
Gaudin Lax matrix and twist function. In addition to the algebra of observables A G 0 , fixed by the realisation π G 0 and the Takiff datum described in the previous subsection, the model we aim to construct depends on the position of the unique site of the AGM, which should be a real number. We shall denote this position by z 1 . The Gaudin Lax matrix of the model is then given by following the general expression (2.15). Similarly, from Equation (2.16), we get the twist function of the model: The zeros of this twist function are If we suppose and ∞ to be of the same sign, these zeros are real and simple (recall that and ∞ are supposed to be non-zero). For future convenience, we will write these zeros as Let us also note that the derivatives of ϕ at these zeros are given by: (3.17) Hamiltonian. From the expressions (3.14) and (3.15) of Γ(z, x) and ϕ(z), one can compute the spectral parameter dependent charge Q(z), defined as (2.35). One then extracts the charges Q 1 and Q 2 associated with the zeros ζ 1 and ζ 2 , using the general definition (2.36). One finds As expected from Subsection 2.4.1, one has Q 1 +Q 2 = P G 0 , with P G 0 the momentum (3.5) of the model.
Let us now discuss the Hamiltonian of the model. Following Subsection 2.2.3, we shall parametrise this Hamiltonian as Moreover, according to Subsection 2.4.2, we shall choose the numbers 1 and 2 to be either 1 or −1 to ensure the relativistic invariance of the model. As discussed above, the choice 1 = 2 = ±1 corresponds to taking H = ±P G 0 : this choice hence leads to a trivial dynamics that we shall not consider. Thus, we should choose either 1 = − 2 = 1 or 1 = − 2 = −1. These two choices lead to opposite Hamiltonians and thus equivalent dynamics. We shall focus on the choice which makes the Hamiltonian positive (if g 0 is a compact real form), which is obtained by taking i of the same sign as −ϕ (ζ i ). By inspecting Equation (3.17), we thus choose 1 = 1 and 2 = −1. The corresponding Hamiltonian is

The model in Lagrangian formulation
Inverse Legendre transform on the fields. Recall that the algebra of observables A G 0 is described by the G 0 -valued field g and the g 0 -valued field X, encoding a field on T * G 0 equipped with its canonical symplectic structure. Thus, the model described above is an Hamiltonian field theory on a canonical cotangent bundle T * G 0 and hence should admit a Lagrangian description as a two-dimensional field theory on G 0 . This Lagrangian description is the subject of this subsection.
To obtain it, one needs to perform an inverse Legendre transform. In particular, one needs to pass from the Hamiltonian fields g(x) and X(x) to a Lagrangian field g(x, t). This is done by expressing the conjugate momenta of the model, encoded in X, in terms of the time derivative of g. For that, let us compute the time evolution of g defined from the Hamiltonian H by

From the expression (3.18) of H and the Poisson brackets (3.3), (3.4) and (3.8), one gets
This gives the Lagrangian expression of the Hamiltonian field X.
Action. Recall that the canonical Hamiltonian field on T * G 0 is encoded in terms of coordinates fields φ i and their conjugate momenta π i . The action of the model is then obtained from the Hamiltonian as the inverse Legendre transform where one should replace all Hamiltonian fields by their Lagrangian expressions. From the definition (3.2) of X, one can rewrite the first term of this action in terms of currents, giving 20) or again One recognizes in the last term of this equation the Wess-Zumino term I WZ [g] as expressed in (3.11). The first term is easily expressed in terms of the Lagrangian fields using Equation (3.19). Similarly, using this equation and the fact that j = g −1 ∂ x g by definition, one obtains the Lagrangian expression of the Hamiltonian (3.18): Combining all these, one gets the expression of the action of the model: We recover here the action of the Principal Chiral Model with arbitrary Wess-Zumino term, as announced. This justifies a posteriori the name PCM+WZ realisation for the Takiff realisation defined in Subsection 3.1.3 (as well as the name PCM realisation in the case where = 0, as the model then coincides with the PCM itself, without Wess-Zumino term).
Lax pair. By construction, the model that we are considering is integrable and in particular, possesses a Lax pair. The general expression of this Lax pair in light-cone coordinates is given by Equation (2.60) in terms of the currents Γ(ζ 1 , x) and Γ(ζ 2 , x). From Equations (3.14) and (3.17), one can express these currents and thus the Lax pair in terms of the currents j and X − W . Using the Lagrangian expression (3.19) of the latter, one then obtains the Lax pair of the model in its Lagrangian formulation. More concretely, one gets the standard light-cone Maurer-Cartan currents.

Free parameters of the model
Redundancy in the parameters. The model initially introduced possessed 4 parameters: the levels , and ∞ and the position z 1 . One easily sees from the Hamiltonian (3.18) or the action (3.21) that there is a redundancy in these parameters. Indeed, the only parameters of the models which appear in these expressions are and K.
Change of spectral parameter. As explained in Subsection 2.3.2, different realisations of AGM can correspond to the same model if they are related by a change of spectral parameter. More precisely, such a change of spectral parameter transforms the levels and the positions of the sites without changing the Hamiltonian. This is the origin of the redundancy in the parameters , , ∞ and z 1 .
Let us be more precise. As explained in Subsection 2.3.2, a change of spectral parameter is given by the combination τ a,b of a dilation by a ∈ R * and a translation by b ∈ R. There are several ways to fix the redundancy of the parameters. One can for example adjust a and b to fix the position of one site and the constant term in the twist function to particular values. In the present case, we will choose to fix the values of the zeros of the twist function to +1 and −1. Thus, we need to find a and b such that This is easily done from the expression (3.16) of ζ 1 and ζ 2 . One gets In what follows we describe the characteristics of the model with spectral parameter z = τ a,b (z).
Twist function and Lax pair. The transformation of the twist function under a change of spectral parameter is given by Equation (2.43). In the present case, we get This twist function indeed has zeros at z equal to +1 and −1. It coincides with the twist function of [20], with the parameter A of [20] equal to 0 (as we are not considering a Yang-Baxter deformation) and with the parameter k of [20] equal to . This is coherent, as the action (3.21) coincides with the one considered in [20] with the Yang-Baxter deformation turned off. Similarly, one can find the Lax pair of the model after the change of spectral parameter using Equation (2.47). Explicitly, this gives

The model in Hamiltonian formulation.
Coupling several PCM+WZ σ-models. In the previous subsection, we identified the PCM with Wess-Zumino term as a realisation of AGM with one site. In Subsection 2.3.3, we developed a natural method to couple realisations of AGM while preserving integrability. Thus, one can apply this method to construct an integrable model coupling N copies of the PCM with Wess-Zumino term.
The description of this model, both in the Hamiltonian and Lagrangian formulations, is the subject of the present subsection. As explained in the previous subsection, the model with a single copy possesses one real site of multiplicity 2, realised in terms of canonical fields in T * G 0 through the PCM+WZ realisation. Following the construction of Subsection 2.3.3, the model with N copies then possesses N real sites of multiplicity 2, that we shall denote (1) to (N ), each associated with a copy of the PCM+WZ realisation. These realisations are parametrised by levels (r) 0 ∈ R and (r) 1 ∈ R * . Similarly to the case N = 1 discussed in the previous subsection, we will use the following simpler notation: Following Subsection 2.3.3, the algebra of observables of the coupled model with N copies is given by A = A ⊗N G 0 . The r th -tensor factor A G 0 in this algebra is described by a canonical field in T * G 0 , as in Subsection 3.1. As explained in this subsection 3.1, this field can be decomposed into a G 0 -valued field g (r) (x) and a g 0 -valued current X (r) (x). Each site (r) of the coupled model is then associated with two real Takiff currents J (r) [0] (x) and J (r) [1] (x), realised as J (r) and J (r) with j (r) = g (r) −1 ∂ x g (r) and W (r) the equivalent of the currents j = g −1 ∂ x g and W of Subsection 3.1 seen in the r th -tensor factor A G 0 in A. The fields g (r) (x), X (r) (x), j (r) (x) and W (r) (x) satisfy the same brackets (3.3), (3.4), (3.8) and (3.9) as the ones in a unique copy of the algebra A G 0 . Moreover, two such fields Poisson commute if they are attached to different sites r. In particular, we shall need in the next subsection the following brackets: For simplicity, we denote by z r = z (r) ∈ R the position of the site (r) in the coupled model with N copies. As explained in Subsection 2.3.3 and more precisely in Theorem 2.1, the decoupling limit of the r th -copy corresponds to taking z r → ∞ while keeping the other positions and all levels fixed.
Gaudin Lax matrix and twist function. Applying the general definition (2.15) of the Gaudin Lax matrix to the present case, one gets Similarly, the twist function of the coupled model is given by It possesses M = 2N zeros ζ i , i ∈ {1, · · · , M }. Finding the expression of these zeros in terms of the levels and the positions of the sites is in general quite complicated, if not impossible, as it requires solving a polynomial equation of order 2N . We will come back to this at the end of this subsection. For now, we will not try to express the ζ i 's explicitly and will just see them as implicitly depending on the levels and the positions through the relation ϕ(ζ i ) = 0. Note that the twist function can be expressed as (3.28) We have shown in Subsection 3.2 that the twist function of the model with one copy has simple and real zeros, if and ∞ have the same sign. For simplicity, we will suppose that for the model with N copies, all levels r , r ∈ {1, · · · , N }, as well as the constant term ∞ , are positive numbers. This ensures, that the twist function of the r th -copy alone possesses simple and real zeros, for all r ∈ {1, · · · , N }. According to Theorem 2.1 and Appendix C, this implies that, at least close enough to the decoupling limit, the zeros of the coupled twist function are also real and simple. More concretely, it means that there is a domain in the space of parameters (with the positions z r large enough), in which these zeros are real and simple. In this domain, we can construct an integrable relativistic coupled model following the construction of Section 2 (which requires the zeros of the twist function to be simple). Moreover, if G 0 is compact, one can choose the parameters of the model such that its Hamiltonian is positive, as all the zeros are real.
Hamiltonian. The Hamiltonian of the model is constructed as explained in Section 2, as the linear combination (2.38) of the quadratic charges Q i with coefficients i equal to ±1 to ensure the relativistic invariance of the model. For the rest of this subsection, we will suppose that there are as many i 's equal to +1 that there are i 's equal to −1 (this is possible as there are M = 2N such i 's) 11 . This will be motivated later, in Subsection 3.3.2, by observing that different numbers of i 's equal to +1 and −1 lead to a model on which we cannot perform an inverse Legendre transform and thus which does not possess a naive Lagrangian formulation. As a consequence of this hypothesis, note that the sets I ± defined in (2.57) are both of size N . This separation of the zeros allows us to write the twist function in a factorised form: . (3.29) The charges Q i are extracted from the spectral dependent quantity Q(z), which is defined in Equation (2.35). Using this definition together with the expressions (3.26) and (3.27) of the Gaudin Lax matrix Γ(z, x) and the twist function ϕ(z) for the model with N copies, one can get an explicit formula for Q(z) and thus for H in terms of the currents j (r) (x) and X (r) (x) − r W (r) (x). One can then compute the time evolution of different observables using ∂ t = {H, ·} and describe the dynamics of the model. This is not the approach that we will develop in this article. Indeed, as we aim to describe integrable σ-models, we will focus more on the Lagrangian formulation, which is treated in the next subsection.

The model in Lagrangian formulation
Inverse Legendre transform on the fields. In the previous Subsection, we described an integrable Hamiltonian model with observables A = A ⊗N G 0 . Similarly to the case with one copy, this model should admit a Lagrangian formulation as a two-dimensional field theory on G N 0 , with fundamental fields g (r) (x, t), for r ∈ {1, · · · , N }. This formulation is obtained from the Hamiltonian one by performing an inverse Legendre transform.
Usually, the first step in this transform would be to express the conjugate momenta of the theory, encoded in the Hamiltonian fields X (r) (x), in terms of the time derivative of the Lagrangian fields g (r) (x, t), through the dynamics ∂ t = {H, ·} of the model. This would require to have an explicit expression of the Hamiltonian H, which we did not seek for in the previous subsection. Instead, we will use a slightly different but simpler approach: as we will see, it will be easier to describe the time derivatives ∂ t g (r) (x, t) in terms of the currents Γ(ζ i , x), i ∈ {1, · · · , M }, instead of the currents X (r) (x). This relation will contain implicitly the inverse Legendre transform on the conjugate momenta and will turn out to be enough to obtain the Lagrangian formulation of the model.
The dynamics of the model is determined by the Hamiltonian H, which itself is expressed in terms of the quadratic charges Q i . In what follows, we shall use the explicit expression (2.37) of the charges Q i in terms of the currents Γ(ζ i , x). To compute the Poisson bracket of g (r) (x) with Q i , let us then first compute its bracket with Γ(ζ i , x). From the expression (3.26) of the Gaudin Lax matrix Γ(z, x) and the Poisson brackets (3.25), one finds g (r) Combining this with the expression (2.37) of Q i , one finds Thus, the Hamiltonian expression of the field j (r) where we used the expression (2.38) of the Hamiltonian. Note that this equation contains the relation between the Lagrangian fields j using the expression (2.52) of the momentum P A . Comparing to Equation (2.39), we then get ± (x, t) = L ± (z r , x, t). Yet, we know from Equations (2.60) that L ± (z, x, t) has simple poles exactly at the ζ i 's, for i ∈ I ± . Moreover, recall from Paragraph 3.3.1 that we supposed that the subsets I ± of {1, · · · , M = 2N } are both of size N . Thus, L ± (z, x, t), seen as a rational function of z, is the sum of N simple fractions. It is a classical result, made precise by Lemma B.2, that such a function is entirely determined by its evaluation at N distinct points. Thus, Equation (3.33) entirely determines L ± (z, x, t) by interpolation. By applying Lemma B.2, one then gets the Lagrangian expression of the Lax pair: where Note that a crucial hypothesis for the reasoning above to work is that the light-cone components L ± (z) of the Lax pair have exactly N simple poles, so that we can interpolate them through their values at the positions z r , r ∈ {1, · · · , N }. This is true as we supposed that there are as many i 's equal to +1 than to −1. Let us suppose that this hypothesis is false: one of the two sets I ± is then of size strictly less than N and the corresponding Lax pair component L ± (z) thus is the sum of less than N simple fractions. It is impossible that the evaluation of such an object at the N points z r 's yields the N linearly independent quantities j (r) ± , hence a contradiction. This observation is an indirect proof that one cannot perform the inverse Legendre transform if one does not have |I + | = |I − | = N (in the sense that the relation between the time derivatives of the Lagrangian fields g (r) (x, t) and the Hamiltonian conjugate momenta cannot be inverted).
Back to the inverse Legendre transform. Recall that the Lax matrix is defined from the Gaudin Lax matrix and the twist function by Equation (2.19). From their expressions (3.26) and (3.27) in the present case, one checks that the current X (r) can be extracted from the Lax matrix through where L (z) denotes the derivatives of L(z) with respect to the spectral parameter z. Using the expression (3.34) of the light-cone Lax pair L ± (z), we find . (3.37) In order to get a presentation as uniform as possible, we will make the choice to express X (r) in terms of ∞ , the z s 's and the ζ i 's only, and thus not to keep explicit dependences on the levels s and s .
To determine the expression of these levels in terms of the z s 's and the ζ i 's, we will use r = χ r (z r ) and r = − together with the formula χ r (z) = − ∞ ϕ +,r (z)ϕ −,r (z), which follows from Equations (3.28) and (3.35). We then get r = − ∞ ϕ +,r (z r )ϕ −,r (z r ) and r = ∞ 2 ϕ +,r (z r )ϕ −,r (z r ) + ϕ +,r (z r )ϕ −,r (z r ) . Finally, we obtain the Lagrangian expression of X (r) , or more precisely X (r) − r W (r) : with ρ rr = ∞ 4 ϕ +,r (z r )ϕ −,r (z r ) − ϕ +,r (z r )ϕ −,r (z r ) , (3.40a) Lagrangian expression of the currents Γ(ζ i ). Let us fix i ∈ I ± . From Equation (2.60), we note that Yet, the Lagrangian expression of L ± (z) was found above using the interpolation lemma B.2. The last statement of this lemma also allows the computation of the residue of L ± (z) at ζ i . More precisely, we get Lagrangian expression of the Hamiltonian. To obtain the Lagrangian expression of the Hamiltonian, let us note that Using the Lagrangian expression (3.42) of the currents Γ(ζ i ), we obtain: .

Thus, we can write
. For r = s, this coefficient can again be rewritten as Thus, we obtain or again For r = s, the coefficient (3.44) becomes Using the identity (3.45), we get or again Action of the model. The inverse Legendre transform is performed at the level of the action as in the case with one copy. Similarly to Equation (3.20) for one copy, we have here for N copies: which can be rewritten as S g (1) , · · · , g (N ) = N r=1 1 2 R×D dt dx κ X (r) − r W (r) , j (r) (3.48) using the expression (3.11) of the Wess-Zumino term. Combining the Lagrangian expressions (3.39) of X (r) − r W (r) and (3.43) of H, we get S g (1) , · · · , g (N ) = 1 2 N r,s=1 R×D dt dx ρ sr κ j (s) where µ ± rs = ρ rs + ρ sr − 4c ± rs 4 .

Parameters of the model.
Positions and levels.
Let us discuss what are the parameters of the model described above. Following the construction of the model in the previous subsections, we see that it is defined by the following parameters: • the positions z r , r ∈ {1, · · · , N }, of the sites ; • the levels r and r , r ∈ {1, · · · , N }, of the sites ; • the constant term ∞ ; • the choice of parameters i = ±1, i ∈ {1, · · · , M }, such that |I There exists a redundancy in these parameters, coming from changes of spectral parameter, as explained in Subsection 2.3.2. More precisely, such changes of spectral parameter remove two degrees of freedom to this set of parameters, without changing its Hamiltonian and action. For example, one can rescale the spectral parameter to fix the constant term ∞ (or one of the levels r ) to a fixed value. Similarly, one can shift the spectral parameter to fix one of the positions z r . In addition to this change of spectral parameter, there exist some discrete redundancies in the parameters: any permutation of the positions and of the corresponding levels leads to an equivalent model (through the same permutation on the fields g (r) ). The choice of parameters described above is particularly adapted to discuss the decoupling limits of the model. Indeed, following Subsection 2.3.3 and Theorem 2.1, one can make the r th -field g (r) decouple from the others by sending the position z r to infinity, while keeping all levels and all other positions fixed.
However this choice of parameters also has downsides. A rather harmless one is the fact that the levels and the positions should be such that the twist function of the model has simple zeros, so that one can apply the construction of the Hamiltonian. Moreover, it is sometimes useful (in particular to ensure the boundedness of the Hamiltonian from below) that these zeros are also real. This is true only on a certain domain of the parameters (see end of Subsection 3.3.1), whose explicit description is quite involved.
The main downside of the present choice of parameters is that the action of the model depends on them only implicitly through the zeros ζ i of the twist function. Yet, the explicit computation of these zeros is in general a very complicated, if not impossible task, as it requires the resolution of a polynomial equation of degree 2N .
Positions and zeros. The observation above shows that the choice of the positions and levels as the parameters of the model is not a very appropriate one. Instead, it suggests to take the positions and the zeros of the twist function as the defining parameters of the model. More precisely, the model is determined by the following parameters • the N positions z r , r ∈ {1, · · · , N }, of the sites ; • the N zeros ζ i , i ∈ I + , whose associated i is equal to +1 ; • the N zeros ζ i , i ∈ I − , whose associated i is equal to −1 ; • the constant term ∞ . This choice of parameters indeed entirely determine the coefficients ρ rs and r appearing in the action (3.49), through Equations (3.38) and (3.40).
The parameters described above are not all free, as they are still subject to the redundancy from changing the spectral parameter. A dilation τ a,0 of the spectral parameter acts on the positions and the zeros by multiplication by a (z r → az r and ζ i → aζ i ) and on the constant term ∞ by multiplication by a −1 ( ∞ → a −1 ∞ ). It is clear from Equations (3.38) and (3.40) that the coefficients ρ rs and r appearing in the action are invariant under such a change. One can for example use this freedom to fix ∞ . A translation τ 0,b of the spectral parameter only affects the positions and the zeros (z r → z r + b and ζ i → ζ i + b) but not the constant term ∞ . For this case also, it is clear from Equations (3.38) and (3.40) that the coefficients ρ rs and r are invariant under such a change.
Note also that there exist discrete redundancies between the parameters. Indeed, any permutation of the zeros {ζ i , i ∈ I + } or the zeros {ζ i , i ∈ I − } does not change the model (permutations mixing the two families of zeros however do yield a distinct model).
This choice of parameters allows a simpler and more direct description of the model than the previous one, as the action of the model is explicitly expressed in terms of these parameters. Moreover, this choice of parameters is also natural to describe the integrable structure of the model. Indeed, the points ζ i define the poles of the Lax pair of the model. Moreover, the Hamiltonian integrability of the model is governed by its twist function, which is easily expressed in terms of the points z r and ζ i as (3.28).
At first, the main downside of this choice of parameters seems to be the discussion of the decoupling limits. Indeed, these limits require to let a position z r flow to infinity, while keeping all levels s and s fixed. This, however, does not correspond to keeping the ζ i 's fixed (in fact, as explained in Appendix C, some of the ζ i 's also go to infinity in this limit). Thus, such a decoupling limit in this parametrisation is difficult to discuss (although, we know that it exists).
However, there exists a slightly different notion of decoupling limit which is more adapted to the parametrisation in terms of positions and zeros. Let us write the number of sites of the model as N = N 1 + N 2 , with N 1 and N 2 positive integers. We introduce a coupling constant γ in the parameters of the model as follows, at first following the initial idea described in Subsection 2.3.3. For r ∈ {1, · · · , N 1 }, we suppose that the positions z r are independent of γ and we write them as z (1) r . For r ∈ {1, · · · , N 2 }, we suppose that the positions z N 1 +r are given by z (2) r + γ −1 . In the "standard" decoupling limit proposed in subsection 2.3.3, the dependence of the zeros ζ i , i ∈ I ± , on the coupling constant γ is such that the levels r and r are fixed. This is where the alternative decoupling limit differs from the standard one. Let us divide the sets I ± labelling the zeros into I ± = I (1) ± I (2) ± , such that |I (1) ± | = N 1 and |I (2) ± | = N 2 (which is always possible as |I ± | = N = N 1 + N 2 ). We will suppose that the zeros ζ i of the model depend on γ as follows: i 's and ζ (2) i 's independent of γ. Then, one checks explicitly that in the decoupling limit γ → 0, one has where the actions S k , for k ∈ {1, 2}, are the actions of the coupled model with N k copies, positions z (k) r (r ∈ {1, · · · , N k }) and zeros ζ ± , associated with the same parameters i = ±1). The decoupling limit at the level of the Lax pair L ± (z) is obtained in a similar way than for the standard one (see Paragraph 2.3.3), by considering the g×g-valued Lax pair L ± (z), L ± (z + γ −1 ) . Note that the introduction of the coupling parameter γ in the parameters of the model can be interpreted nicely in terms of the factorisation of the twist functions and as one then has

Symmetries of the model
Conserved charge associated with a PCM+WZ realisation. For this paragraph and the next one, let us come back to a more general realisation of local AGM, as described in Section 2. We suppose that there exists a fixed real site α 0 ∈ Σ r which is associated with a PCM+WZ realisation, as described in Subsection 3.1 in terms of canonical fields g(x) and X(x) in the cotangent bundle T * G 0 . Thus, the Takiff currents of the site α 0 are given by with the fields j(x) and W (x) defined as in Subsection 3.1. Let α ∈ Σ \ {α 0 } be a site different from α 0 and p ∈ {0, · · · , m α − 1}. The Takiff current J α [p] Poisson commutes with the currents J α 0 [0] and J α 0 [1] . In particular, it is in involution with j and hence with g. We thus deduce that it is also in involution with W and finally with X, using the above expression for J α 0 [0] . As a conclusion, we have: (3.50) Let us define the following g 0 -valued current and charge: One shows that the Poisson brackets of the current K α 0 (x) with the Takiff currents J α 0 [0] (x) and J α 0 [1] (x) are simply: Integrating these brackets over y ∈ D, we obtain As the current K α 0 also Poisson commutes with the Takiff currents associated with other sites α ∈ Σ \ {α 0 } by Equation (3.50), we get In particular, the charge K α 0 is in involution with the quadratic charges Q i defined by (2.37), or more generally the charges Q d i of the integrable hierarchy (2.41). As the Hamiltonian of the model is a linear combination of the Q i 's, the charge K α 0 is conserved: Global symmetry associated with a PCM+WZ realisation. One checks that the current K α 0 (x) is a Kac-Moody current of level 2 : Thus, the conserved charge K α 0 satisfies the Kirillov-Kostant bracket associated with the Lie algebra g 0 : K α 0 1 , K α 0 2 = C 12 , K α 0 1 . By the Noether theorem, this conserved charge is associated with an infinitesimal g 0 -symmetry of the model, defined on any observable O ∈ A by with ∈ g 0 the infinitesimal parameter of the transformation. From Equation (3.51), we see that the Gaudin Lax matrix is invariant under this symmetry: From Poisson brackets (3.3), (3.4) and (3.8), one checks that the symmetry δ α 0 corresponds to the left multiplication of the field g(x): This infinitesimal transformation lifts to the global G 0 -symmetry g → h −1 g of the model, with global parameter h ∈ G 0 .
Global symmetries of the integrable σ-model with N coupled copies. Let us now come back to the integrable coupled σ-model described in this subsection. It is constructed from N independent copies of the PCM+WZ realisation. According to the two previous paragraphs, each of this realisation is associated with a global G 0 -symmetry of the model. Thus, the model possesses a G N 0 -symmetry. The corresponding transformation of the fields g (r) (x, t) of the model is simply with h 1 , · · · , h 1 ∈ G N 0 the constant parameters of the symmetry. We will refer to this as the left G N 0 -symmetry of the model. The global charges associated with this symmetry are simply given (in the Hamiltonian formulation) by Using the Lagrangian expression (3.39) of the field X (r) − r W (r) , as well as the expression (3.38) of r , one finds that the Lagrangian expression of the conserved charge K (r) . More precisely, one finds that its density K (r) (x) coincides with the temporal component 1 2 (K (r) It is obvious that the light-cone currents j The model also possesses another natural symmetry, which comes from the diagonal symmetry present in all realisations of AGM, as discussed in Subsection 2.2.4. It is easy to check that in the present case, the global conserved charge associated with this symmetry is The corresponding symmetry acts by right multiplication on all the fields g (r) (x, t): with h ∈ G 0 . We will refer to this symmetry as the diagonal right G 0 -symmetry of the model. It acts by conjugacy on the light-cone currents: The invariance of the action (3.49) is then ensured by the ad-invariance of the bilinear form κ. It is interesting to note that, due to the non trivial coupling ρ rs between different copies r = s, the right multiplication of the fields g (r) (x, t) by different parameters h r ∈ G 0 is not a symmetry of the action. This symmetry is however recovered in the decoupled limit, as expected from the symmetries of the PCM plus Wess-Zumino term alone.
Using the Lagrangian expression (3.39) of the field X (r) − r W (r) one can write the Lagrangian expression of the density K ∞ (x) of the conserved charge K ∞ . After a few manipulations, one finds that it is equal to the temporal component 1 2

Other realisations and deformations
In this Section, we present other suitable Takiff realisations than the PCM+WZ one. Except for the non-abelian T-dual realisation, these realisations share the same observables A G 0 as the PCM+WZ one. One can use these realisations to construct other relativistic integrable field theories as realisations of local AGM. The simplest examples of such theories correspond to taking these realisations alone.
As we shall see and as was shown in [4], these theories can be identified with other known integrable σmodels: the homogeneous and inhomogeneous Yang-Baxter deformations of the PCM, the non-abelian T-dual of the PCM and its λ-deformation.
As in the case of the PCM+WZ realisation, one of the interest of understanding these models as realisations of local AGM is that one can then couple them, using the general method developed in Subsection 2.3.3 of this article. This results in a very rich landscape of deformations of the integrable coupled σ-model constructed in Subsection 3.3. We do not aim to describe and analyse this entire landscape here, as this would require many further investigations. Instead, we present the general ideas allowing the construction of such deformed models and give a simple example to illustrate how this general method works.

The realisation
Let us consider the algebra of observables A G 0 , parametrised by the canonical fields g(x) and X(x), respectively in G 0 and g 0 , as described in Subsection 3.1. Let us also consider a solution R : g 0 → g 0 of the homogeneous Classical Yang-Baxter Equation (CYBE) on g 0 : which we suppose to be skew-symmetric with respect to the non-degenerate form κ: For ∈ R * , we define the following g 0 -valued currents [3]: The skew-symmetry (4.2) of R implies that Using this equation and the CYBE (4.1) and starting with the Poisson brackets (3.3) and (3.4), one checks that Thus, J [0] (x) and J [1] (x) are real Takiff currents with multiplicity 2 and levels 0 = 0 and 1 = .
They define a Takiff realisation in A G 0 , which we call the homogenous Yang-Baxter (hYB) realisation. This realisation has the same multiplicity and levels that the PCM realisation (without the current W (x), i.e. for = 0, see Subsection 3.1.3). In fact, it is easy to check from Equation (3.12) with = 0 that these two realisations coincide for R = 0. Thus, the hYB realisation is a deformation of the PCM one through the additional introduction of a solution R of the CYBE.
The generalised Segal-Sugawara integrals of the hYB realisation can be obtained using the formula (2.9) giving their expression for a general Takiff realisation of multiplicity 2. One then gets where we used the skew-symmetry (4.2) of R and recognized the expression (3.5) of the momentum P G 0 of the algebra A G 0 . Thus, the hYB realisation is suitable, as defined in Paragraph 2.1.4.

The model with one copy of the realisation
As we defined a suitable Takiff realisation with observables A G 0 , one can construct an integrable relativistic field theory on A G 0 as a local AGM in this realisation. The corresponding AGM will then be composed of only one site of multiplicity two. As the levels of the hYB realisation are the same as the one of a PCM realisation, the twist function of the model coincides with the one (3.15) with = 0. As explained in Subsection 3.2.3, there is a redundancy in the parameters of the twist function (3.15), coming from changes of the spectral parameter. One way to fix this redundancy was to fix the values of the zeros of the twist function to +1 and −1. For simplicity, we shall make the same choice for the model with the hYB realisation. The twist function is then given by similarly to (3.22) (with = 0). The corresponding Gaudin Lax matrix is then Similarly to the PCM case, we attach to the zero ζ 1 = +1 the parameter 1 = +1 and to the zero ζ 2 = −1 the parameter 2 = −1. One can then compute the Hamiltonian (2.38) of the model, by finding the expression of the quadratic charges Q 1 and Q 2 from the twist function and the Gaudin Lax matrix above. We get From this Hamiltonian, one computes the dynamics of the field g(x) and finds the Hamiltonian expression of j 0 = g −1 ∂ t g. Inverting this expression, we get in terms of the light-cone currents j ± = g −1 ∂ ± g. One can then perform the inverse Legendre transform and find the action of the model, yielding One recognizes in this expression the action of the homogeneous Yang-Baxter deformation of the PCM, first introduced for the AdS 5 × S 5 superstring in [21]. Note that the diagonal symmetry present for all realisations of local AGM corresponds here to the right multiplication on the field g(x, t).

Homogeneous Yang-Baxter deformations
Let us consider a realisation of local AGM, with arbitrary sites and Takiff realisations. Let us suppose moreover that one of the real sites α 0 ∈ Σ r is associated with a PCM realisation (without current W (x), i.e. with α 0 0 = 0). One can then construct another integrable model, defined as a realisation of local AGM with exactly the same twist function, but with the site α 0 realised in terms of the hYB realisation instead of the PCM one. By construction, this model is an integrable deformation of the one initially considered, which now depends on the matrix R (it is a deformation in the sense that the limit R = 0 gives back the initial model). This is the principle of homogeneous Yang-Baxter (hYB) deformations. Of course, the simplest example of such a deformation is the one constructed in the previous paragraph: it is the hYB deformation of the PCM itself, seen as a realisation of local AGM with only one site.
However, this method now allows the construction of a wider class of hYB deformed models. In particular, one can apply this deformation to the integrable coupled σ-model with N copies described in Subsection 3.3, as long as it possesses a site (r) with corresponding level (r) 0 = −2 r = 0. One can even apply a hYB deformation to all sites with r = 0, yielding a multi-deformed model involving different R-matrices. As an example, we develop in Appendix D the case where we apply a deformation to only one site of the model.

The realisation
Let R : g 0 → g 0 be a linear map on g 0 , skew-symmetric with respect to the bilinear form κ, i.e. satisfying Equation (4.2). We also suppose that it is a solution of the modified (inhomogeneous) Classical Yang Baxter Equation (mCYBE): with c = 1 (we then talk of a split R-matrix) or c = i (we then talk of a non-split R-matrix).
Let us consider the algebra of observables A G 0 , described by canonical fields g(x) and X(x) on T * G 0 , as described in Subsection 3.1. Forγ a non-zero real number, we define the following currents [2,3] J (±) with R g = Ad −1 g • R • Ad g , as in (4.3). If we are considering the split case c = 1, then the currents J (±) [0] (x) are real, in the sense that they are g 0 -valued. If we are considering the non-split case c = i, then the currents J [0] 2 (y) = 0, using the mCYBE (4.7) and the skew-symmetry (4.4) of R g . Thus, the currents J (±) [0] (x) are two commuting Kac-Moody currents (i.e. Takiff currents of multiplicity 1), with levels .
According to the discussion above, they are real Kac-Moody currents in the split case and complex ones in the non-split case.
As the currents J (±) [0] (x) are Kac-Moody currents satisfying the appropriate reality conditions, they define a Takiff realisation of the Takiff algebra T with Takiff datum defined as follows. If c = 1, the Takiff datum has two real sites (+) and (−) in Σ r , with multiplicity 1 and levels (±) 0 = ±(2γ) −1 , and no complex sites (Σ c = ∅). If c = i, the Takiff datum has no real sites (Σ r = ∅) and one complex site (+), with multiplicity 1 and level As these sites are of multiplicity one, they are associated with a unique generalised Segal-Sugawara integral, whose general expression is given in Equation (2.8). In the present case, using the skewsymmetry (4.2) of R, one easily checks that with P G 0 the momentum (3.5) of A G 0 . Thus the Takiff realisation defined above is suitable. We call it the inhomogeneous Yang-Baxter (iYB) realisation.

The model with one copy of the realisation
The iYB realisation defined above can be used to construct relativistic integrable models using the general method of Section 2. As for the PCM+WZ realisation (see Subsection 3.2) and the hYB realisation (see Subsection 4.1), let us first describe the model obtained when considering a realisation of local AGM with just one copy of the iYB realisation.
This model possesses two sites (+) and (−), either real in the split case, or complex conjugate in the non-split one. These sites are associated with positions z (+) and z (−) , either real or complex conjugate. Recall that one has the freedom to shift all the positions of the sites by the same real number without changing the model (via a translation of the spectral parameter). One can use this freedom to require z (+) = −z (−) : we will denote by cη the value of z (+) (we then have η real in both the split and non-split case). The twist function of the model is therefore given by This twist function possesses two opposite zeros. We shall suppose that they are real. By a dilation of the spectral parameter, one can then always send these zeros to +1 and −1, so that the deformed model shares the same zeros as the PCM (see Equation (3.22) with = 0). This, up to a redefinition of the parameter η, amounts to fix ∞ to an appropriate value. In the end we get The Gaudin Lax matrix of the model is then given by (taking into account thatγ = η K ) As in the undeformed case (or the homogeneous Yang-Baxter deformation), we attach to the zero ζ 1 = +1 the parameter 1 = +1 and to the zero ζ 2 = −1 the parameter 2 = −1, thus fixing the Hamiltonian (2.38) of the model. From the expressions of the twist function and Gaudin Lax matrix above, one compute the quadratic charges Q 1 and Q 2 using Equation (2.37) and thus the Hamiltonian. One gets Note the similarity with the Hamiltonian (4.6) of the homogeneous Yang-Baxter deformation. One can then perform the inverse Legendre transform on this Hamiltonian and compute the action of the deformed model. One finds One recognizes the action of the inhomogeneous Yang-Baxter deformation of the PCM [22,23]. As for the hYB deformation of the PCM, the global diagonal symmetry is associated with the right multiplication on the field g(x, t).

Inhomogeneous Yang-Baxter deformations
Construction. As in the homogeneous case discussed in Subsection 4.1, let us consider a realisation of local AGM which possesses a real site α 0 ∈ Σ r attached to a PCM realisation. The twist function and Gaudin Lax matrix of the model are thus of the form where ϕ(z) and Γ(z, x) contain the informations about the parameter ∞ and the other sites α ∈ Σ \ {α 0 } of the model. Let us consider the "deformed twist function" It has simple poles at z α 0 ± cη, with residues res z=zα 0 ±cη It is thus the twist function of a local AGM with sites obtained from the initial model by replacing the site α 0 (of multiplicity 2, levels α 0 1 and α 0 0 = 0 and position z α 0 ) by two sites α ± 0 (of multiplicity 1, levels α ± 0 0 = ± 1 2cγ and positions z α ± 0 = z α 0 ± cη) 12 . We will suppose that the sites α ∈ Σ \ {α 0 } keep the same realisations as in the initial model. Moreover, we shall realise the two sites α ± 0 through the inhomogeneous Yang-Baxter realisation described above. The "deformed" Gaudin Lax matrix is then: One checks that the deformed twist function ϕ η (z) and deformed Gaudin Lax matrix Γ η (z, x) gives back the undeformed one ϕ(z) and Γ(z, x) in the limit η → 0. Thus, we have constructed a deformation of the initial model with deformation parameter η. By construction, this deformation is integrable for all values of η. This is the principle of what we call inhomogeneous Yang-Baxter deformations, which thus provides a general scheme to construct integrable deformation of local AGM which possess a PCM realisation. The simplest example of such a deformation is the deformation of the PCM itself, seen as a realisation of local AGM with a single site. Up to a redefinition of the parameters, this simplest deformation coincides with the one (4.8) presented in the previous paragraph.
The deformation method presented here only changes the characteristics of the model which are associated with the site α 0 (it does not modify the levels and positions of the other sites or the constant term ∞ ). Although this seems like the most natural way of introducing the deformation, it also has some downsides. The main one is that the modification (4.9) of the twist function changes its zeros. In particular, the zeros of ϕ η (z) cannot be expressed in terms of the zeros of ϕ(z) in a easy way. As these zeros are one of the key ingredients in the construction of the integrable deformed model, this can make the description of the deformation quite involved in general.
One can already see this in the inhomogeneous Yang-Baxter deformation of the PCM alone: the parametrisation that we choose to describe this model in a simple way in the previous paragraph does not correspond exactly to the one we would get from the general deformation scheme introduced here. Indeed, in the previous paragraph, we choose to deform the twist function in a way which preserves its zeros, which then also requires to modify appropriately the constant term ∞ in the twist function. This model is however equivalent to the one we would obtain from the general deformation scheme (which would let ∞ unchanged but which would modify the zeros of the twist function), through a redefinition of the parameters.
This illustrates that there is not a unique way to apply inhomogeneous Yang-Baxter deformations. The general scheme presented here is the only one which can be applied to all local AGM with a PCM realisation, regardless of the other sites of the model. However, in some cases, it seems easier to introduce a deformation of the twist function which allows a simpler description of its zeros (by ensuring that they are not modified or at least that they can be easily expressed). Whether all these deformations are equivalent up to a redefinition of the parameters seems at the moment like a difficult question and would require further investigations.
It is natural to apply inhomogeneous Yang-Baxter deformations to the integrable coupled σ-model of Subsection 3.3. In principle, one can construct many multi-parameter deformations of this model, as one can apply the general deformation scheme to different copies of the model (as long as the corresponding level (r) 0 = −2 r is equal to 0). The description of the whole landscape of these deformations is out of the reach of the present article, as we just aim to describe the general idea at the origin of the construction.
Inhomogeneous Yang-Baxter deformations as q-deformations. It is a well-known fact [2] (see also [24,25,26]) that the inhomogeneous Yang-Baxter deformation of the PCM breaks the left multiplication symmetry of the model and changes it into a q-deformed symmetry (when we are considering a model on the real line D = R). We shall now show that this a general feature of the inhomogeneous Yang-Baxter deformation scheme introduced in the previous paragraph. We shall use the notations of this paragraph: in particular, the undeformed model possesses a PCM realisation, attached to the site α 0 . Recall from Subsection 3.3.4 that this implies the existence of a global symmetry in the undeformed model, corresponding to a left-multiplication g → h −1 g of the field g ∈ G 0 attached to the site α 0 , which generalises the left symmetry of the PCM alone. It is easy to check that this symmetry is broken by the deformation scheme introduced above.
Let us now consider the Lax pair L η (z, x), M η (z, x) of the deformed model. One can apply to it a formal gauge transformation, which preserves the zero curvature equation. In particular, let us apply such a gauge transformation with gauge parameter the field g itself: We will be more precisely interested in the evaluation of this gauge-transformed Lax matrix at the poles z α ± 0 of the twist function (which are the positions of the sites α ± 0 of the deformed model): One can evaluate the field L η z α ± 0 , x using Equation (2.22). We get . with It is explained in [2,3,27,7] how Equation (4.10) implies the existence of a q-deformed Poisson-Hopf algebra of conserved charges of the model, with q = e −icγ , extracted from the monodromy matrices of the currents L g η z α ± 0 , x (which are conserved if D = R). Moreover, these q-deformed conserved charges are associated with a Poisson-Lie symmetry of the model, whose explicit expression can be found in [27,7] and which gives back the global left symmetry g → h −1 g in the limit η → 0.

Inhomogeneous Yang-Baxter deformations with Wess-Zumino term
In the previous subsections, we presented iYB deformations of PCM realisations, without Wess-Zumino term. However, it is known [20] that there also exists an inhomogeneous Yang-Baxter deformation of the PCM with Wess-Zumino term (at least in the non-split case and when R is the standard Drinfel'd-Jimbo R-matrix). According to the section 3.4 of the article [20], this deformed model is characterised by two commuting Kac-Moody currents (obtained from the evaluation of the Lax matrix at the poles of the twist function). These currents are a generalisation of the ones J (±) [0] introduced in this subsection, with an additional dependence on the current W (x). The main difference with the currents J (±) [0] is that the Kac-Moody currents of [20] do not possess opposite levels. They define a generalisation of the iYB realisation described in this subsection.
Similarly to the deformation scheme described above, this more general realisation allows the deformation of local AGM which possess a PCM+WZ realisation. In particular, this generalised deformation scheme can be applied to all copies of the integrable coupled σ-model of Subsection 3.3, without any assumptions on the undeformed model. For brevity and simplicity, we choose not to develop this more general deformation here.

The realisation
In this subsection, we will consider the phase space A g 0 generated by a coordinate field v(x) valued in g 0 and its conjugate momentum. The coordinate field is described by scalar fields v a (x) by expanding it in the basis {I a } a∈{1,··· ,n} of g 0 : v(x) = v a (x)I a .
Each of the fields v a (x) possesses a conjugate momentum m a (x). We encode these conjugate momenta in the following g 0 -valued field m(x) = κ ab m a (x)I b .
The canonical Poisson brackets For ∈ R * , let us define the following g 0 -valued currents From the Poisson brackets (4.11), we get Thus, the currents J [0] (x) and J [1] (x) are Takiff currents of multiplicity 2 and levels 0 = 0 and 1 = and hence define a Takiff realisation in A g 0 , that we shall call non-abelian T-dual realisation, for reasons to be made clear below.
The generalised Segal-Sugawara integrals D [0] and D [1] of this Takiff realisation can be computed from the general formula (2.9). In particular, we get One checks that this quantity generates the spatial derivatives ∂ x on the fields v(x) and m(x) and thus coincides with the momentum P g 0 of the algebra A g 0 . The non-abelian T-dual realisation is hence suitable.

The model with one copy of the realisation
One can construct an integrable relativistic model as a realisation of a local AGM with one site, associated with the non-abelian T-dual realisation described above. Via a change of spectral parameter, the twist function of this model can be chosen to be as for the PCM with no Wess-Zumino term (see Subsection 3.2.3). The corresponding Gaudin Lax matrix is As for the PCM, we associate the zeros ζ 1 = +1 and ζ 2 = −1 of the twist function ϕ(z) with the parameters 1 = +1 and 2 = −1. One can then compute the Hamiltonian (2.38) of the model, finding To perform the inverse Legendre transform of this model, we first express the momenta m in terms of the time derivative ∂ t v = {H, v}. Rewriting the result in terms of the light-cone derivatives Note in particular that this implies (4.14) The action of the model is then given by the inverse Legendre transform where we replace m and K ∂ x v + [m, v] by their Lagrangian expression (4.13) and (4.14). Using the skew-symmetry of ad v with respect to the bilinear form κ, one finds One recognizes here the action of the non-abelian T-dual of the PCM [28,29]. Let us note that the global diagonal symmetry discussed in §2.2.4 acts here on v(

Equivalence with the PCM realisation
The non-abelian T-dual model obtained above is equivalent to the PCM. More precisely, it was shown in [30] that the two models are related by a canonical transformation, i.e. a Poisson-preserving isomorphism from the algebra A G 0 to the algebra A g 0 . Let us consider a realisation of local AGM which possesses a site α 0 associated with a PCM realisation. One can consider another realisation of AGM with the same sites and levels, but with the PCM realisation replaced by the non-abelian T-dual one. These two models are then related by a canonical transformation (which changes the algebra A G 0 attached to the site α 0 into the algebra A g 0 and leaves the other Takiff realisations unchanged) and are thus equivalent. This construction can be applied in particular to the integrable coupled σ-model of Subsection 3.3, applying this transformation to part or to all of the N copies of the algebra A G 0 (under the additional condition that the corresponding levels r are equal to 0). For the sake of brevity, we will not develop this model further here.

The realisation
Let us consider the algebra A G 0 of canonical fields in T * G 0 , as described in Subsection 3.1. Let us consider [31,3,32,33] the g 0 -valued current (4.16) and its "gauge transformation" [0] (x) considering the canonical fields (g −1 , −gXg −1 ) instead of (g, X) and taking to − . Note that it coincides with the density of the conserved charge associated with a PCM+WZ realisation, as introduced in Paragraph 3.3.4 (we already observed in that paragraph that it is a Kac-Moody current and that it Poisson commutes with the Kac-Moody current of the PCM+WZ realisation, i.e. with J (+) [0] (x)). The two Kac-Moody currents J (±) [0] (x) together define a Takiff realisation with two sites of multiplicity one, in the algebra A G 0 , which we call the λ-realisation.
The generalised Segal-Sugawara integrals D (±) [0] of the λ-realisation can be computed from the general expression (2.8). In particular, using the identity (3.10), one finds with P G 0 the momentum (3.5) of the algebra A G 0 . Thus, the λ-realisation is suitable.

The model with one copy of the realisation
As for the other realisations above, let us describe the local AGM which one obtains if we consider only one copy of the λ-realisation. The twist function of this model possesses two simple poles, whose residues are opposite one to another (as the λ-realisation is obtained from two Kac-Moody currents with opposite levels). Via a change of spectral parameter, one can always choose this twist function to have zeros at +1 and −1, similarly to all the other cases described above. The twist function can then be parametrised as for some non-zero real numbers α and K (the notation ϕ λ (z) is introduced here for further convenience, as we shall use later a new parameter λ instead of α). This function can be written as This partial fraction decomposition allows us to read the positions z ± = ±α and the levels (±) 0 = ∓2 of the two sites (+) and (−) of the model. The Gaudin Lax matrix of the model is then given by with J (±) [0] (x) the Kac-Moody currents (4.16) and (4.17) of the λ-realisation. From the expressions above of the twist function and the Gaudin Lax pair of the model, one can compute the quadratic charges Q 1 and Q 2 associated with the zeros ζ 1 = +1 and ζ 2 = −1. The Hamiltonian H of the model is given from these zeros by the linear combination (2.38), where we choose 1 = +1 and 2 = −1 (as in the other cases of this section). After a few manipulations, one finds where we introduced λ = 1 + α 1 − α .
Note that, in terms of λ, the level reads From the expression of the Hamiltonian, one can compute the time evolution of the field g(x) and more precisely the current j 0 = g −1 ∂ t g in terms of the currents j and X − W . One can then invert this relation and find the Lagrangian expression of the Hamiltonian field X − W . Here, we give this expression in terms of the light-cone currents j ± = g −1 ∂ ± g: One can then perform the inverse Legendre transform of the model and compute its action. In the end, one finds where S WZW, [g] is the Wess-Zumino-Witten action We recognize in the expression (4.18) the action of the λ-deformed model [34]. The symmetry g(x, t) → h −1 g(x, t)h of this action corresponds to the global diagonal symmetry present for all local AGM.

λ-deformations
Let us consider a realisation of local AGM which possesses a real site α 0 ∈ Σ r attached to a nonabelian T-dual realisation. The twist function and Gaudin Lax matrix of the model are thus of the form where ϕ(z) and Γ(z, x) contain the informations about the parameter ∞ and the other sites α ∈ Σ \ {α 0 } of the model. Let us consider the "deformed twist function" It has simple poles at z α 0 ± α, with residues res z=zα 0 ±α It is thus the twist function of a local AGM with sites obtained from the initial model by replacing the site α 0 (of multiplicity 2, levels α 0 1 and α 0 0 = 0 and position z α 0 ) by two sites α ± 0 (of multiplicity 1, levels α ± 0 0 = ∓2 and positions z α ± 0 = z α 0 ± α). We will suppose that the sites α ∈ Σ \ {α 0 } keep the same realisations as in the initial model. Moreover, we shall realise the two sites α ± 0 through the λ-realisation described above. The "deformed" Gaudin Lax matrix is then: It is easy to see that, in the limit α → 0, one has Thus, one can expect that the model with twist function ϕ α (z) is a deformation (with deformation parameter α) of the model initially considered, with twist function ϕ(z). To prove this statement, one would need a similar limit with the corresponding Gaudin Lax matrix Γ α (z, x) of the model. However, it is easy to see that the naive limit α → 0 cannot be taken in the Gaudin Lax matrix (4.21). This is in contrast to the analogous limit in the context of the inhomogeneous Yang-Baxter deformation. This limit requires therefore a more subtle treatment. We will suppose [34] that in the limit α → 0 we have where v(x) is thus a g 0 -valued field. One then checks that We shall also suppose that: for some g 0 -valued field m(x). One can then read the Poisson brackets of the fields v(x) and m(x) from the ones (3.3) of the fields g and X, by considering the limit α → 0. We find that the fields v(x) and m(x) satisfy the same bracket (4.11) as the fields of the same introduced in Paragraph 4.3.1 and that appears in the undeformed Gaudin Lax matrix (4.19). The limit α → 0 of the deformed Gaudin Lax matrix (4.21) can be computed from the asymptotic expansions (4.22), (4.23) and (4.24), noticing that In the ends one finds as expected. Thus, the construction of this paragraph allows to define, for any local AGM which possesses a non-abelian T-dual representation, an integrable deformation of this model with deformation parameter α. The simplest example of such a deformation is the deformation of non-abelian T-dual model alone, which is the model presented in the previous paragraph.

Conclusion
The reformulation of the class of integrable field theories admitting a twist function as affine Gaudin models was proposed in [4] as a natural framework for addressing the problem of quantisation of such models. In this article we have shown that this language can already be put to good use at the classical level to build new classical integrable field theories with infinitely many parameters. One of the virtues of this approach is that, by working directly at the Hamiltonian level, the classical integrability of the resulting field theory is automatic. A usual drawback of working at the Hamiltonian level is that Lorentz invariance is not explicit, making it more difficult to identify which of the integrable field theories constructed are relativistic. However, it turns out that there are quite simple conditions ensuring Lorentz invariance for the class of models we have considered.
It is worth emphasising that one of the key properties of Gaudin models, which lies at the heart of the construction presented in this article, is the fact that they are 'ultralocal' as spin chains, in the sense that observables associated with different sites mutually Poisson commute. Indeed, in the affine Gaudin model description of an integrable field theory the fundamental fields are combined into particular Takiff currents, thought of as attached to individual sites, which commute with one another. This is exemplified in both the Yang-Baxter deformation and λ-deformation of the principal chiral model, both of which are realisations of affine Gaudin models with two sites to which are attached mutually Poisson commuting Kac-Moody currents.
This property of 'ultralocality' is the one which ultimately makes it possible to assemble together integrable field theories by treating them as affine Gaudin models. We gave a detailed description in Subsection 3.2 and Section 4 of various integrable field theories that can be used as building blocks in such a construction. We illustrated the method for assembling together an arbitrary number of these building blocks on two examples, in Subsection 3.3 and Appendix D. In turn, the resulting theories can themselves be used as building blocks in subsequent constructions.
Since our starting point is the affine Gaudin model description of an integrable field theory, the construction also makes apparent how the parameters of the affine Gaudin model, which are encoded in the twist function, enter in the corresponding action and the Lagrangian expressions for the Lax pairs of these models.
It is clearly very desirable to apply the method presented in this paper to assemble together other integrable σ-models in a variety of different ways, in order to determine the full scope of possibilities. It will then also be interesting to explore the infinite family of integrable field theories arising from this construction. We expect that some integrable σ-models which have appeared recently in the literature belong to this infinite family. This should, for instance, be the case for the coupled integrable λ-models constructed in [35,36] (see also [37]). It is expected that the integrable field theory constructed even more recently in [38] is also an affine Gaudin model.
On another front, it would be interesting to study the renormalisation group flow for the integrable σ-models obtained by using the present construction. A rich structure has already been found in the renormalisation group flow of the integrable λ-models considered in [35,36].
To borrow terminology used in the context of Gaudin models, the ones considered in this article are all 'non-cyclotomic'. A natural follow-up will therefore be to generalise the present construction to the family of integrable field theories which are realisations of cyclotomic affine Gaudin models. This family includes [4] the symmetric space σ-models and, more generally, integrable Z T -coset models [39]. In fact, more importantly, all these theories have only a regular singularity at infinity, which is responsible for their gauge symmetry. Specifically, the absence of an irregular singularity at infinity gives rise to a first-class constraint [4] which generates the gauge invariance in the model. Another example of a field theory where this happens is a gauge-invariant version of the bi-Yang-Baxter model which was described in [40] as a one-parameter deformation of the η-deformed Z 2 -permutation coset space σ-model. As an affine Gaudin model it is non-cyclotomic and has only a regular singularity at infinity [4].
Recall, by contrast, that the integrable field theories considered in this article all have a double pole at infinity. More precisely, the double pole manifests itself only in the form of a constant term − ∞ in the partial fraction decomposition of the twist function ϕ(z). The effect of such a term, which is clear from (2.42), is to break the gauge symmetry inherent in affine Gaudin models with a regular singularity at infinity down to the global symmetry discussed in Paragraph 2.2.4.
It will be interesting to understand how the method of assembling integrable σ-models presented here can be generalised to the case when the individual theories have gauge symmetries. In particular, since first-class constraints are always associated with the same point, namely infinity, there will still only be one single constraint in the assembled theory which will now correspond to a diagonal gauge symmetry. Specifically, as already announced in [1], one may expect that starting from a Lie group G 0 and a subgroup H of G 0 such that the σ-model on the coset G 0 /H is integrable, it is possible to construct a G ×N 0 /H diag σ-model which is integrable. This certainly requires a detailed investigation.
A.2 Extracting the generalised Segal-Sugawara integrals from Q(z) Let us consider the charge Q(z) as defined in (2.35). The following Proposition shows that one can naturally extract the generalised Segal-Sugawara integrals D α [p] = π D α [p] (in the realisation π : T → A) from Q(z). Proposition A.3. One has where Q reg is regular at all positions z α .
Proof. Let us fix a site α ∈ Σ. Using Equations (2.15) and Lemma A.2, we can write: From the definition (2.35) of Q(z), we then get In this expression, we have −2m α + 2 ≤ q + r − s ≤ m α − 1. The terms with −2m α + 2 ≤ q + r − s ≤ −1 are regular at z = z α and can thus be reabsorbed in O (z − z α ) 0 . We can then write where we extended the sum over 0 ≤ s ≤ 2m α − 2 as η α s = 0 for s ∈ {0, · · · , m α − 1}. For fixed p, q and r, the index s is fixed to the value q + r − p by the presence of the term δ p,q+r−s , under the condition that 0 ≤ q + r − p ≤ 2m α − 2. As, we have p ≥ 0 and q, r ≤ m α − 1, the second inequality in this condition is always fulfilled. Thus, we get This ends the demonstration, as one gets by applying the realisation π to Equation (2.6).
We now have enough to prove Proposition A.1.
Proof. From the definition (2.26) and (2.35) of H(z) and Q(z), we have By Proposition A.3 and the definition (2.23) of P(z), we then get This proves Proposition A.1, as ϕ(z) has a pole of order m α at z = z α and Q reg (z) is regular at z = z α .

B Interpolating rational functions
Lemma B.1. Let V be a vector space, let N and m 1 , · · · , m N be positive integers and let z 1 , · · · , z N be pairwise distinct complex numbers. Let us consider the following V -valued rational function of a complex parameter z: where the v r,p 's, for r ∈ {1, · · · , N } and p ∈ {0, · · · , m r −1}, are elements of V . Let M = m 1 +· · ·+m N and ζ 1 , · · · , ζ M be pairwise distinct complex numbers, also pairwise distinct with the z r 's. Then the v r,p 's are linear combinations of the f (ζ i )'s.
Proof. Let us define the following polynomial It is clear that g(z) = α(z)f (z) is a polynomial of degree at most M − 1. It is a classical result that such a polynomial is uniquely determined by its evaluation at the M distinct points ζ i . More precisely, we have with λ i (z) the Lagrange interpolation polynomials. We then have For r ∈ {1, · · · , N } and p ∈ {0, · · · , m r − 1}, we have hence the lemma.
We will also need the following lemma, which is a more precise version of the one above in the case where all poles of the function f (z) are simple. Note however that we change some of the notations between the two lemmas (in particular we exchange the roles of the z r 's and the ζ i 's), to adapt to the different contexts in which they are used in the main text of the article.
Lemma B.2. Let V be a vector space, let N be a positive integer and let z 1 , · · · , z N , ζ 1 , · · · , ζ N be 2N pairwise distinct complex numbers. Let us consider the following V -valued rational function of a complex parameter z: where the v i 's, for i ∈ {1, · · · , N }, are elements of V . Then the function f (z) is uniquely determined by its evaluation at the N points z r , for r ∈ {1, · · · , N }. More precisely, one has the following interpolating formula: .
It is also clear that h(z) is a polynomial of degree at most N − 1. As it vanishes at N different points, it is identically zero, which proves We then get the last statement of the lemma by where we used ϕ r (ζ i ) = 0. for sites α both in Σ 1 and Σ 2 . These currents form the realisation π 1⊗2 of the Takiff algebra T 1⊗2 . The way we define the level at infinity ∞ and the positions w α of the sites of the coupled model is explained in Paragraph 2.3.3. Following this definition, the twist functions of the decoupled models are

C Coupling and decoupling of realisations of local AGM
for k equal to 1 or 2, while the twist function of the coupled model is (C.1) Similarly, the Gaudin Lax matrix of the coupled model is given by Decoupling limit. Let us consider the decoupling limit γ → 0, as defined in Subsection 2.3.3. It is easy to check that Thus, the twist function and Gaudin Lax matrix of the coupled model go to the ones of the first model in the decoupling limit. This asymmetry between the two decoupled models can seem surprising at first. It actually comes from the way we introduced the coupling parameter γ in Equation (2.50), as a shift of the positions of the second model with respect to the positions of the second one, which we keep fixed. One can also consider the same model by keeping the position of the second model fixed while shifting the positions of the first. This is done by translating the spectral parameter by γ −1 . For this new spectral parameter, the twist function and Gaudin Lax matrix of the coupled model go to the ones of the first model in the decoupling limit, i.e.

C.2 Zeros of the coupled twist function
Decoupling limit of the zeros. In this paragraph, we will need the following lemma, which is a classical result of algebraic geometry: Lemma C.1. Let f (γ) (z) = 0 be an algebraic equation on a complex number z, depending algebraically on a complex parameter γ. If x(0) is an isolated solution when γ = 0, then for sufficiently small γ, there exists an isolated solution x(γ), depending algebraically on γ, such that x(γ) = x(0) + O(γ).
We then have f 1⊗2 ζ i (γ) = 0. Thus, ζ i (γ) is a simple zero of the coupled twist function ϕ 1⊗2,γ (z), such that Considering the asymptotic expansions (C.4) and (C.5) and the fact that the zeros ζ (k) i are simple, it is clear that the zeros ζ i (γ) of ϕ 1⊗2,γ (z) are all different, at least for γ small enough. As there are M of them, they exhaust the whole set of zeros of ϕ 1⊗2,γ (z). In conclusion, we have proved the following result, which makes the first part of Theorem 2.1.
but differs through its Gaudin Lax matrix, which is now given by with R (N ) = R g (N ) = Ad −1 g (N ) • R • Ad g (N ) . The Hamiltonian H hYB of the deformed model is defined as usual as the linear combination (2.38) of the quadratic charges Q i . In order to construct an actual deformation (which gives back the initial model when R = 0), we choose the parameters i of the deformed Hamiltonian to be the same as the undeformed one. As in the undeformed case, we shall not try to describe explicitly the Hamiltonian of the model and instead will focus on finding its Lagrangian formulation, starting with the expression of the Lax pair.

D.3 Lagrangian Lax pair of the model
Lagrangian Lax pair from interpolation. Recall that in the undeformed case, we found the expression of the Lagrangian Lax pair by interpolation, using mainly the fact that the evaluation of the Lax pair L ± (z) at the position z = z r yields the light-cone current j (r) ± (see Equation (3.33)). Going through the proof of this result in Subsection 3.3.2, one checks that it just relies on the facts that the site (r) is associated with the PCM+WZ realisation. Thus, this proof also applies in the present deformed case for the sites (1) to (N − 1), which are still associated with a PCM+WZ realisation, hence: L hYB ± (z r ) = j (r) ± , ∀ r ∈ {1, · · · , N − 1}.
This property suggests to use a similar interpolation method as the one we used in the undeformed case. Recall that this method could be applied because the light-cone components L ± (z) of the Lax pair possess exactly N simple poles (the zeros ζ i for i ∈ I ± of the twist function). This property is still true for the deformed case as we did not change the coefficients i . Thus, similarly to Equation 1 (x), Γ hYB (ζ i , y) = From the expression (2.37) of Q hYB i and the skew-symmetry of R (N ) , we get Thus, the Hamiltonian expression of the field j (N ) 0 = g (N )−1 ∂ t g (N ) is given by where M hYB (z, x) denotes the derivative of M hYB (z, x) with respect to the spectral parameter z. A similar computation with all i 's replaced by 1, i.e. with H replaced by the momentum P A , yields the same relation with j (N ) (x) and L hYB (z N , x). Thus, we also have j (N ) This is the generalisation for the deformed case of the relation (3.33).
Finding the current J

D.4 Inverse Legendre transformation and action of the model
Inverse Legendre transform of the fields. To perform the inverse Legendre transform of the model, one needs to find the relation between the conjugated momenta, encoded in the currents X (r) , and the time derivatives of the coordinates fields, encoded in the currents j (r) 0 = g (r)−1 ∂ t g (r) . In the undeformed case, this was done in Equation (3.39), using the identity (3.36). The proof of this identity relies solely on the fact that the site (r) is attached to a PCM+WZ realisation. Thus, it is also true in the present deformed case for the sites (1) to (N − 1). This allows one to extract the Lagrangian expression of X (r) − r W (r) , for r ∈ {1, · · · , N − 1}, in a very similar way than Equation (3.39) in the undeformed case. The main difference with the undeformed case is that the derivatives L ± (z r ) appearing in the identity (3.36) have now to be computed with the expression (D.5) of the deformed Lax pair L hYB ± (z), where the current j i.e. that Equation (D.11) is also satisfied for r = N , as N = 0.
Lagrangian expression of the currents Γ hYB (ζ i ) and the Hamiltonian. Following the method developed in the undeformed case, we know turn to the Lagrangian expression of the current Γ hYB (ζ i ). One easily checks that the computations of the undeformed case stay true in the present one if one replaces the currents j with the coefficients c ± rs given by Equations (3.46) and (3.47).
Action in terms of J (r) ± 's. Let us finally turn to the computation of the action. The expression (3.48) of the inverse Legendre transform is also valid for the deformed case (remembering that now N = 0). Thus, one can write S g (1) , · · · , g (N ) = Inserting this in the expression for the action above, we get S g (1) , · · · , g (N ) = R×D dt dx for all r, s ∈ {1, · · · , N }. Furthermore, we note that for one of the indices equal to N , this expression simply reduces to O rN = ρ rN 1 − ϑ R (N ) and O N r = ρ N r 1 − ϑ R (N ) , for r ∈ {1, · · · , N }.
Let us conclude by pointing out that it follows directly from Equation (D.5) that the field equations of the deformed model are obtained from those of the undeformed model after replacing the current j