New integrable boundary conditions for the Ablowitz-Ladik model: from Hamiltonian formalism to nonlinear mirror image method

Using Sklyanin's classical theory of integrable boundary conditions, we use the Hamiltonian approach to derive new integrable boundary conditions for the Ablowitz-Ladik model on the finite and half infinite lattice. In the case of half infinite lattice, the special and new emphasis of this paper is to connect directly the Hamiltonian approach, based on the classical $r$-matrix, with the zero curvature representation and B\"acklund transformation approach that allows one to implement a nonlinear mirror image method and construct explicit solutions. It is shown that for our boundary conditions, which generalise (discrete) Robin boundary conditions, a nontrivial extension of the known mirror image method to what we call {\it time-dependent boundary conditions} is needed. A careful discussion of this extension is given and is facilitated by introducing the notion of intrinsic and extrinsic picture for describing boundary conditions. This gives the specific link between Sklyanin's reflection matrices and B\"acklund transformations combined with folding, {\it in the case of non-diagonal reflection matrices}. All our results reproduce the known Robin boundary conditions setup as a special case: the diagonal case. Explicit formulas for constructing multisoliton solutions on the half-lattice with our time-dependent boundary conditions are given and some examples are plotted.


Introduction
The problem of formulating integrable partial differential equations (PDEs) on finite or half-infinite intervals appeared very soon after the discovery of the Inverse Scattering Method (ISM) [1,2], in [3]. However, this did not address, in the context of boundary conditions, one of the most important aspects of integrable systems: the fact that integrable PDEs (in 1 + 1 dimensions) are (infinite-dimensional) Hamiltonian systems. Sklyanin's seminal work [4] paved the way to a framework that tackles boundary conditions in integrable classical systems both from the point of view of ISM and of Hamiltonian theory 1 . The PDE point of view seems to have prevailed however until recently and developed into the so-called nonlinear mirror image method, following the initial impetus of [5,6,7,8] and revived more recently e.g. in [9,10,11]. Note also a new related angle on the question, called boundary dressing, which appeared in [12,13]. The key tools of that method are Bäcklund transformation and a special folding symmetry.
On the one hand, the Hamiltonian approach to integrable boundary conditions rests upon the classical reflection equation which describes the allowed boundary conditions via its solutions: the reflection matrices. The latter are conditioned by the classical r-matrix [14,15] of the model under consideration. On the other hand, the nonlinear mirror image approach relies on the use of a Bäcklund matrix which satisfies what we could call a boundary zero curvature equation. Both aspects are contained in [4]. However, a thorough analysis of the class of solutions allowed by one or the other method shows that the reflection equation admit more general solutions (non-diagonal reflection matrices) than the boundary zero curvature equation. This begged the question: why such a discrepancy and how can we resolve it? This question was the motivation behind [16] where a detailed answer was provided. The non-diagonal reflection matrices can be cast into a boundary zero curvature form provided a time dependent term is included. This time-dependent term is in fact the time derivative of the reflection matrix itself. In other words, the way out of this discrepancy is to consider a time-dependent version of reflection matrices. This was of course well-known in the Hamiltonian setup where they appear as solution of the so-called dynamical reflection equation [17]. As it turns out, the time-dependent boundary zero curvature equation has also appeared before in the literature but its direct connection with Hamiltonian setup was not understood. Rather, some a posteriori checks were performed in some instances (see e.g. [18] and references therein for this point of view).
The main aim of the present paper is to follow the general results obtained in [16] and illustrate them in full detail on the example of the Ablowitz-Ladik (AL) model [19]. In doing so, we actually find new integrable boundary conditions for the AL model. In fact, we work with a more general Ablowitz-Ladiktype model which depends on three arbitrary parameters, see e.g. [20]. For special values of the parameters, one can recover the usual AL model or a discrete modified Korteweg-de Vries model. Our approach is as follows. Starting from the Hamiltonian approach, we obtain the general non-diagonal solution of Sklyanin's classical reflection equation and show how to derive the (boundary) zero curvature equations. Our results contain the case of (discrete) Robin boundary conditions treated previously in [21] as a particular case. Using the method introduced in [22], we compute the Lax pair for the model on the finite lattice. We then make the explicit connection with the nonlinear mirror image method for our new integrable boundary conditions. To facilitate the transition, we introduce the notion of intrinsic and extrinsic picture for describing integrable boundary conditions. They are related by a local gauge transformation which realises the transition from the time-independent form of the boundary zero curvature equation (intrinsic picture) to its time-dependent form (extrinsic picture). Once in the extrinsic picture, we are able to use the ideas of the nonlinear mirror image method and to implement them in our case. Explicit formulas for multisolitons solutions in the presence of our new integrable boundary conditions are provided via symmetries on the scattering data coming from the folding procedure.
In Section 2, we review the classical r-matrix approach for the Ablowitz-Ladik (AL) model with periodic boundary conditions in order to introduce the required notations and tools. Note that we obtain different results from the standard ones as we use a more convenient normalisation for the Lax matrix of AL. We recall the notion of monodromy and transfer matrices and how one can derive the Lax pair and the zero curvature representation of the equations of motions from the Hamiltonian formalism. In Section 3, the necessary modifications to the r-matrix approach, as proposed in [4], are implemented for AL on the finite chain with integrable boundary conditions in our normalisation. We derive a general reflection matrix, obtain the corresponding equations of motion on the open interval, first as Hamiltonian equations of motion and then as a zero curvature equation for an appropriately modified Lax pair, to take the boundary conditions into account. In section 4, we introduce the notion of intrinsic and extrinsic picture for boundary conditions. It is first illustrate for the case of (discrete) Robin boundary conditions and then implemented for our more general case. The time-dependent form of our new boundary conditions is then obtained. In Section 5, the Lax pair and (boundary) zero curvature representation of the model with boundary conditions are derived from the Hamiltonian approach. The extrinsic picture is described in this setup, which allows us to make the connection with the nonlinear mirror image method. In Section 6, the results of the previous section are used to construct explicit solutions, in particular multisoliton solutions (see Proposition 6.3), for the problem on the half-lattice with our integrable boundary conditions. We restrict our attention to the discrete NLS reduction in the focusing regime and implement the nonlinear mirror image method. Some conclusions are gathered in the last section.

Transfer matrix formalism and integrable Hamiltonian
The classical r-matrix approach provides a neat and powerful formalism to present the Hamiltonian formulation of integrable systems. It also allows one to derive the time Lax matrix and the zero curvature representation as Hamilton's equations. In this section, we follow the well-known method (see e.g. [23]) and adapt it to the Ablowitz-Ladik model. We stress that the explicit results are actually new since we use a different normalisation of the Lax pair, which turns out to change the r-matrix underlying the model. In particular, the Liouville integrability of the general AL system (2.22)-(2.23) below is established for the first time here to the best of our knowledge. The normalisation we use was studied in [20] where it was advocated to be better than the traditional one as it produces a Lax pair with appealing algebraic properties. We will see throughout this paper that indeed, this normalisation is superior for various reasons.
The starting point of the classical r-matrix approach to Hamiltonian integrable lattice models is the so-called Lax matrix (j, z), where j is an integer associated to a site in the chain and z is the spectral parameter, which is assumed to obey the following quadratic ultralocal Poisson algebra which we consider in C 2 ⊗ C 2 here. This Poisson algebra encodes the Poisson brackets used on the phase space of the model. The standard auxiliary space notation has been used here where I is the 2 × 2 identity matrix, E mn is the canonical basis of 2 × 2 matrices and mn (j, w) are the entries of the matrix (j, z). The classical r-matrix r ab is a 4 × 4 matrix. For the Ablowiz-Ladik model we take Notice the extra factor in front of the matrix which ensures that det (j, z) = 1. With this choice, one can check that having the r-matrix in (2.1) of the form is equivalent to the following Poisson brackets on the fields These are the standard AL brackets. The r-matrix (2.4) is skew-symmetric and satisfies the classical Yang-Baxter equation These two properties ensure that the Poisson bracket defined by (2.1) is antisymmetric and satisfies the Jacobi property. In addition to these properties, the r-matrix (2.4) is symmetric in the auxiliary spaces There is a one-parameter family of normalisations for (j, z): for all s ∈ R. It is still possible to find an r-matrix for each s (2.10) The case s = 0 corresponds to the normalisation which is most often used for the AL model. In the present paper, we use s = −1/2 which is the particular value for which r (s) simplifies nicely and enjoys the additional property (2.8) of being symmetric. This is another argument in favour of this particular normalisation, in addition to those advocated in [20].
The Hamiltonian of the model on N + 1 sites (and all conserved quantities) can be obtained by defining the so-called single-row monodromy matrix L a (z) = a (N, z) a (N − 1, z) . . . a (1, z) a (0, z) , (2.11) and the associated single-row transfer matrix Then, one shows that the following holds 14) The second relation allows us to take t(w) as the generating function in w of elements I (n) in involution: Let us introduce Then, by direct calculation, one shows that where we have used the conventions q N +1 ≡ q 0 and r N +1 ≡ r 0 . Let us introduce the Ablowitz-Ladik type Hamiltonian as the following combination of elements in involution We associate to the Hamiltonian H a time evolution according to With this choice, the equations of motion introduced in [20] appear as Hamilton's equations for the fields q j , r j contained in (j, z) ∂ t (j, z) = {H, (j, z)} , (2.21) i.e. explicitly, for j = 0, 1, . . . , N , together with periodic boundary conditions q N +1 ≡ q 0 , r N +1 ≡ r 0 , q −1 ≡ q N and r −1 ≡ r N . Since these equations of motion derive from the Hamiltonian H which was constructed from the transfer matrix t(z), this proves the Liouville integrability of the model with periodic boundary conditions. The Ablowitz-Ladik equations are recovered for α = β = 1 2 and γ = −1. With the additional reduction r j = νq * j (ν = ±1), equations of motion (2.22) and (2.23) becomes the ones of the (integrable) discrete nonlinear Schrödinger equatioṅ For α = −β = i 2 and γ = 0 with the reduction q j = νr j (and ν = ±1 and q j real), we recover the equations of motion of the discrete modified Korteweg-de Vries equation given bẏ

Lax pair associated to the Ablowitz-Ladik chain
It is one of the remarkable features of the r-matrix approach that instead of guessing a Lax pair for a given system of equations, one can derive the time Lax matrix A(j, z) from the knowledge of the space Lax matrix (j, z) and its associated r-matrix 2 . The zero-curvature form of the equations of motion (or compatibility condition of the Lax pair) is also a by product of this approach. We implement this for the AL chain in the present normalisation. In particular, we will derive the Lax matrix A(j, z) given in [20]. Let us define the partial monodromy for n ≥ m L a (n, m, z) = a (n, z) a (n − 1, z) . . . a (m, z) . (2.26) We use the convention L(n − 1, n, z) = I and obviously one gets L(N, 0, z) = L(z). Following [15,24], one defines, for j = 0, . . . , N + 1, To simplify the expansion in terms of the spectral parameter w, we introduce a regularized version of The matrix M (j, w, z) is the generating function in w of the matrices M (n) (j, z) In particular, one gets Using relation (2.1), one proves Comparing with the expansion (2.15) for t, this shows how the matrices M (n) (j, z) are associated to the charges I (n) . Indeed, by expanding relation (2.34) with respect to w, one gets Remark 2. We see that the explicit form of the matrices M (−N −1) (j, z) and M (N +1) (j, z) given by (2.30) and (2.33) are different whereas the corresponding charges I (−N −1) and I (N +1) are equal to C. This apparent contradiction is solved by noting that we can always add to any M (n) (j, z) a matrix proportional to the identity matrix and independent of the site j. In the following, we will use M (−N −1) (j, z) or M (N +1) (j, z) for the charges C.
In particular, in view of the expression (2.19) of H in terms of the I (n) 's, we obtain Using the Leibniz rule and relation (2.35), this becomes where we have defined The last term in (2.38) is irrelevant in (2.37) but allows us to obtain that A(j, z) be traceless. By using the explicit forms of the matrices M (n) (j, z) given by (2.30)-(2.33), we obtain where σ 3 is the Pauli matrix and This concludes our derivation of the time Lax matrix A(j, z) given in [20]. Equation (2.37) now yields the zero curvature representation of (2.21) In other words, the pair ( (j, z), A(j, z)) is a Lax pair for the system under consideration.

Bäcklund transformation
It is well-known that the existence of the Lax pair ( (j, z), A(j, z)) satisfying the zero curvature equation allows one to write a consistent auxiliary problem for the 2 × 2 matrix φ(j, t, z): For later purposes, it is important to notice that two solutions of the equations of motion can be associated by a Bäcklund transformation. Such a transformation has a representation at the level of the Lax pair (gauge transformation) and at the level of the wavefunction (Darboux transformation). For convenience, we will stick to the name Bäcklund transformation for all these different aspects of these transformations. Suppose that we have two solutions q j , r j andq j ,r j such that there exists a matrix B(j, t, z) satisfying B(j, t, z) φ(j, t, z) = φ(j, t, z) , (2.44) where φ(j, t, z) stands for the wavefunction of (2.42)-(2.43) associated with the fieldsq j ,r j . The matrix B(j, t, z) realises a Bäcklund transformation. The consistency of relation (2.44) with the auxiliary problems (2.42)-(2.43) for φ(j, t, z) and φ(j, t, z) leads to the following equations for B(j, t, z) It is easy to see that the system of the above equations is equivalent to equation (2.45) for all j and equation (2.46) only for one given position j 0 , since ( (j, z), A(j, z)) and (˜ (j, z), A(j, z)) satisfy the zero curvature condition. Then, to obtain the admissible Bäcklund transformation, we must first solve equation (2.45).
The following Lemma gives one such solution that will be used below. We will deal with (2.46) in Section 5.2.
where s j = 1−r j q j 1− r j q j and the fields are constrained by Proof. The proof is done by direct computation by inserting expression (2.47) into (2.45) and matching the powers of the spectral parameter z.
In the solution for the Bäcklund matrix given in the previous lemma, there are 8 free parameters which we can take to be for example f 1 0 , f 2 0 , g 1 0 , g 2 0 , x 1 0 , x 2 0 , y 1 0 and y 2 0 . These are determined by fixing B(0, t, z).
and we recover the result of [21], up to slight modifications due to the different choice of normalisation of the matrices (j, z).

.1 Double-row transfer matrix
To consider open finite chains with integrable boundary conditions, one needs to modify the single-row method of the previous section. In his seminal work, Sklyanin [4] proposed to consider the so-called double-row transfer matrix instead of the single-row transfer matrix (2.11). It is defined by where L a (z) is defined as in (2.12). The matrices k ± (u) are reflection matrices describing the boundary conditions at both ends of the finite chain and τ is a function of the spectral parameter that depends on the model. In our case, Let us emphasize that τ is chosen so that ω(τ (z)) = ω(z) where ω is given by (2.40) and that τ (τ (z)) = z. In the historical paper [4], a simpler involution τ (z) = 1/z was used but we must generalize that construction to find integrable boundary conditions for any α and β. The matrices k − (z) and k + (τ (z)) are required to be solutions of the reflection equation These equations ensure that L a (z) satisfies the so-called dynamical reflection equation which in turn allows one to prove the the analog of the important result (2.14) for b(u) i.e.
We take the following solutions of (3.3) as our reflection matrices where a, b, c and d are arbitrary parameters. The matrix k − (z) is the most general solution of the reflection equation (up to an irrelevant overall scaling). As we shall see, it describes the boundary conditions at the "origin". The matrix k + (z) describes the vanishing of the fields at the site N + 1 and is chosen for simplicity here. In particular, this choice allows us to consider the "half-line" problem with vanishing fields at infinity, simply by taking N to infinity. Those choices are already general enough to obtain new integrable boundary conditions for the Ablowitz-Ladik model on the half-line which contain the (discrete) Robin boundary condition as a particular case. The expansion of the double-row transfer matrix b(z) (3.1) is written as follows From relation (3.5), we deduce that {I (n) , I (p) } = 0. Upon inspection of the bulk terms, we define the Hamiltonian by the following combination of the previous charges: By using the explicit form of b(z), we can show that the Hamiltonian is given explicitly by Comparing with (2.19), we see that the effect of the boundary is contained in B and involves the two neighbouring sites j = 0, 1 of the origin. From now on, the time evolution is associated to this Hamiltonian H by Then, the Hamilton's equations of motion can be computed using (2.5). Explicitely, on the left boundary they reaḋ while in the bulk they are given by, for j = 2, 3 . . . , N − 1, and on the right boundary, we havė The equations of motion on the finite lattice come from the Hamiltonian (3.8) which was extracted from the transfer matrix b(z). This gives a proof that they are integrable in the sense of Liouville. We now turn to the derivation of the time Lax matrix producing these on the finite open lattice.

Double-row Lax pair
Recall that using the r-matrix approach, we can derive the time Lax matrices for all the time flows associated to the Lax matrix (j, z) via its single-row transfer matrix, as well as the corresponding zero curvature equations reproducing Hamilton's equations of motion. The same holds true for models with boundaries but, of course, one has to modify formulas (2.27) and (2.28) to take into account the presence of integrable boundaries, as dictated by the double-row transfer matrix. Following the construction of [22], we define, for j = 0, 1, . . . , N + 1, These matrices satisfy the following property [22] Note that the proof of [22] must be adapted here to take in account our involution τ . In [16], it was shown that M(j, w, z) satisfies two other equations which are of crucial importance for describing integrable boundary conditions in a zero curvature representation Writing the expansion of M(j, w, z) as where ω(z) is given by (2.40). The zero curvature equations (3.25)- (3.27) are invariant by adding a term to A(j, z) which is proportional to the identity, is independent of j and is invariant under z → τ (z). We use this freedom in (3.28) to add iω(z) and make A(j, z) traceless. Then ( (j, z), A(j, z)) is an adequate Lax pair associated with the AL model with integrable boundary conditions determined by the matrices k ± (z). The status of relations (3.26)-(3.27) is discussed in detail in [16]. Upon performing the explicit computations of the matrices M (0) (j, z), M (1) (j, z), we find that the matrices A(j, z) for j = 2, . . . , N have exactly the same form as the ones with periodic boundary condition. Namely, one gets However, the matrices A(0, z) and A(1, z) are different and Similarly, the matrix A(N + 1, z) is also different from the bulk one and is given by  The reader familiar with integrable chain models on a finite interval will be perfectly content with the intrinsic representation of the boundary conditions at the end of the interval within the equations of motion as in (3.12)- (3.19). However, the reader who is more familiar with integrable PDEs on the finite interval (or half-line) might be more used to an extrinsic representation of the problem in the form of a bulk equation of motion valid for all values of the space coordinates supplemented by a condition on the field (and spatial derivatives) at the coordinate of the boundary.
To clarify what we mean, let us first consider (3.12)- (3.19) in the case c = d = 0, which corresponding the (discrete) Robin condition on the left boundary as we will see. The intrinsic picture is given by eqs (3.12)- (3.19) which boil down to the equations in the bulk, now valid for j = 1, . . . , N , the equations at the left boundary, and the ones at the right boundary, We see that the equations at j = 0 and at j = N are different from the naive continuation of the bulk equations for j = 1, . . . , N − 1 to j = 0 or j = N . In the extrinsic picture, the idea is to enforce this continuation and to think of the boundaries as sitting at j = −1 and j = N + 1. To produce a system of equation equivalent to the one above, we introduce the value of the fields q −1 , r −1 , q N +1 and r N +1 at these sites such that the bulk equations of motion are the same for all j = 0, . . . , Ṅ Then, at j = 0, the equivalence with the intrinsic picture is restored when imposing the boundary conditions In the discrete NLS case α = β = 1/2, these are known to be the discrete analog of the Robin boundary conditions, as studied e.g. in [21]. At j = N , the equivalence with the intrinsic picture is obtained by setting q N +1 = 0 , r N +1 = 0 , (4.10) which are Dirichlet boundary conditions. In the rest of the paper, we keep these conditions on the right boundary and send N → ∞ to consider the system on the half-infinite lattice with zero boundary conditions at infinity. We also restore c, d = 0 in order to consider our general boundary conditions in the extrinsic picture.
The inversion of the change of variables (4.11) gives In these new variables, the equations of motion for j ≥ 1 reaḋ Then, by using the equations of motion (3.12) and (3.13), the previous relations can be writteṅ In view of this, the equations of motion at the site 0 for Q 0 and R 0 can be written as the equation of motion in the bulk (4.13)-(4.14) continued to j = 0, i.e. : provided Q −1 and R −1 satisfy the following conditions We can summarize the previous discussion in the following proposition.  Let us remark that the boundary conditions (4.21)-(4.22) depends on three sites -1,0 and 1. However, by using equations of motion (4.23)-(4.24) at site 0, we can eliminate Q 1 and R 1 to get boundary conditions depending only on the two sites -1 and 0. The time derivative of the fields at the site 0 then appears. Namely, the boundary conditions (4.21)-(4.22) become Due to the presence of the time derivative in the boundary conditions, we call them time-dependent boundary conditions. This terminology has in fact a deeper root as will become clear later on: the reflection matrix describing them is time-dependent. This has non-trivial consequences on the implementation of the so-called nonlinear mirror image method which we address below.

Reductions
The reduction of the previous equations of motion to get DNLS or DMKdV leads to some constraints on the boundary parameters. For the DNLS (i.e. α = β = 1/2, γ = −1 and R * j = νQ j ), the extrinsic equations of motion becomeṡ with the boundary conditions The parameters of the boundary satisfy a, b ∈ R, c = −νd * . For the DMKdV (i.e. α = −β = i/2, γ = 0, R j = νQ j and R j ∈ R), the extrinsic equations of motion becomesQ with the boundary conditions The parameters a, b, c, d of the boundary should all be real and satisfy c = −νd, b = 0.

15
In view of the previous discussion, we need to establish the effect of going from intrinsic to extrinsic picture on the results of Section 3.2.

Change of variable, gauge transformation and extrinsic picture in the Lax pair presentation
By using the change of variable (4.11), we get new expressions for the matrices A(j, z) for j ≥ 1, We see in particular that the Lax matrix A(1, z) on site j = 1 takes the same form as the bulk ones A(j, z), j ≥ 2. The same feature appears already in Section 4.1 where the equations of motion on the site j = 1 becomes similar to the ones of the bulk after the change of variable. However, even after the change of variables, the matrix A(0, z) still has a different structure from the bulk ones. This is not compatible with the spirit of the extrinsic picture and we have seen that introducing fields Q −1 and R −1 satisfying appropriate boundary conditions (see (4.21)-(4.22)) allowed us to have equations of motion with the same form for all j ≥ 0. Therefore, we would like to introduce the same fields such that the Lax matrix A(0, z) is written as in the bulk, i.e.
and such that the equation is equivalent to the boundary conditions (4.21)-(4.22). However, for the generic boundary conditions associated to k − (z), this procedure fails since A(0, z) given by (3.31) cannot be written in the form (5.2). This is readily seen from the difference in the dependence on the spectral parameter z which cannot be accommodated by a constraint involving only fields. This is a feature of our new boundary conditions. Indeed, for the particular choice of the Robin boundary conditions (c = d = 0), the passage to the extrinsic picture actually works and equation (5.3) with A(0, z) given by (5.2) is equivalent to the Robin conditions (4.9). To overcome this problem for generic boundary condition, we consider a gauge transformation G(z) concentrated at site 0 and defined by with, for j ≥ 0 We find that the gauge transformation is given in terms of the fields q 0 and r 0 by Let us emphasize that since the gauge transformation depends on the fields q 0 and r 0 , its time derivative does not vanish. By injecting the gauge transformation for L (j, z) into the zero curvature equation (3.25), one gets that the Lax matrices A(j, z) must be transformed to The remarkable result is that (after some tedious computations) we can show that A (0, z) can indeed be written as a Lax matrix of the same form as the bulk ones and takes the form (5.2) where Q −1 and R −1 satisfies (4.21) and (4.22), as desired. Therefore, we have succeeded in performing the transition from intrinsic to extrinsic picture consistently at the level of the Lax pair: ( (j, z), A(j, z)) → ( (j, z), A (j, z)) .
It is explicitly given by In turn this has a nontrivial impact on the boundary zero curvature equation (3.26) which now reads This equation is the main reason for calling these boundary conditions time-dependent: (5.11) is the timedependent generalisation of (3.26). Of course, it is not the first time that this equation appears in relation to integrable boundary conditions. However, as explained in [16] and as illustrated on the AL model here, it is the first time that it is directly linked to the Hamiltonian approach and that it appears as a necessity to give an extrinsic picture of boundary conditions corresponding to the most general solutions of the reflection equation. Finally, we can now give the extrinsic form of the equations of motion in the zero curvature representation.
Proposition 5.1. The equations of motion given by relations (4.13)-(4.14) for j ≥ 0 and the boundary conditions (4.21) and (4.22) are equivalent to Proof. The procedure explained previously provides the proof that the equations of motion imply the zero curvature representation. The implication in the other way is proven by direct computation.
In other words, the boundary matrix K − (z) is now dynamical and describes time-dependent boundary conditions, even though the original boundary matrix k − (z) was non-dynamical. We wish to stress at this point that this observation was one of the main motivation behind [16]. We see that if one want to consider the most general solution of the non-dynamical reflection equation (3.3) and have an extrinsic interpretation of the boundary conditions, one is naturally led to consider the dynamical setting. It is at the basis of the long-standing discrepancy between the classical r-matrix approach to boundary conditions and the zero curvature approach. This was resolved in general in [16] and this paper is an explicit illustration of the process on the AL model.
Another outstanding question at this stage is how to deal with our time-dependent boundary conditions in the scheme of the so-called nonlinear mirror image. The latter has been well-known since the early papers [5,6,7,8] where it appeared in connection with the use of Bäcklund transformation together with a folding procedure. However, up to now, this scheme has only been use in connection with the time-independent boundary zero curvature equation of the form (3.26). We develop the time-dependent nonlinear mirror image method in two steps. In the next subsection, we present the first step in implementing this approach with the view of constructing explicit solutions: we embed (5.11) into a Bäcklund transformation scheme coupled with a folding prescription. Then, in Section 6, this is used in conjunction with the standard ISM framework to construct explicit solutions of our AL model on the half-line with integrable boundary conditions.

Auxiliary problem, Bäcklund transformation and nonlinear mirror image method
The zero curvature equations presented in Proposition 5.1 provide the consistency relations for the following auxiliary problem for the 2 × 2 matrices Φ(n, t, z) 3 Φ(j + 1, t, z) = L (j, z) Φ(j, t, z) , for j = 0, 1, . . . (5.14) Here, ρ(z) is a function of the spectral parameter independent of n and t, chosen such that The existence of such function can be shown but in the following we do not need its explicit expression. We now show how this auxiliary problem arises from a Bäcklund transformation by using the nonlinear mirror image method. We start from two solutions (Q j , R j ) and ( Q j , R j ) of the AL model on the full line (j ∈ Z) such that they are related by the following folding conditions It is straightforward to show that this transformation leaves the equations of motion of the AL model invariant. We have also the following properties on the Lax pair due to (5.18) In particular, one gets A (0, z) = A (0, τ (z)) . (5.22) Due to this symmetry on the Lax pair, we deduce that Φ(j, t, z) and J j Φ(−j, t, τ (z)) are eigenvectors of the same auxiliary problem. Therefore, there exists a matrix M (z) independent of the time and of the position such that Φ(j, t, z) = J j Φ(−j, t, τ (z)) M (z) .
This is nothing but relation (5.16) and we conclude that the solution Q j and R j satisfy the desired boundary conditions encoded in (5.13). This shows that the folding procedure coupled with the Bäcklund transformation approach yields solutions of the equations of motion with the desired integrable boundary conditions. It turns out that it also implies the following symmetry relation on B(j, t, z) Lemma 5.1. The Bäcklund matrix B(j, t, z) in the folding procedure explained above satisfies Proof. By definition of B(j, t, z), we have and (5.28) is a consequence of (5.30).

19
The results of Lemma 2.1 are still valid in this case, simply replacing q j , r j by Q j , R j . By using condition (5.26), one fixes the arbitrary parameters as As mentioned previously, we have that det(L (j, z)) = 1 which implies that det(B(j, t, z)) is independent of j (it is easy to see that on relation (2.45)). Therefore, we get in particular det(B(+∞, t, z)) = det(B(0, t, z)) = det(K − (z)) .
Note that these determinants are also time independent due to the tracelessness of A and A . To determine the form of B(+∞, t, z), we assume that the fields vanish at infinity (this could be shown following the same argument as in [21]). Looking at the coefficients z 2 and z −2 in relation (5.33), we deduce that and In this case, we obtain that B(+∞, t, z) is in fact independent of t and we have The coefficients of z 0 in relation (5.33) gives a constraint for f 1 ∞ which reads For any solution f 1 ∞ of (5.37), one gets that B(+∞, τ (z))B(+∞, z) = 1 ρ(z)ρ(τ (z)) I. Then, in particular, one gets ϕ(z)ϕ(τ (z)) = 1 ρ(z)ρ(τ (z)) . We want to stress that in this paper, the choice of the value of B(0, t, z) is dictated by the solution K − of the reflection equation. This is in contrast with previous approaches, e.g. [21], where educated guesses for B(0, t, z) are taken in order to produce the desired boundary conditions under the folding method. Our approach has an advantage in practice (no guess work) but it also has a deeper significance which is the underlying message in this paper: the Hamiltonian approach and the zero curvature approach to integrable boundary conditions are not separated topics but two faces of the same coin. This interplay is well-known for problems on the line but somehow was not as clear in the case with boundaries, despite the seminal work of Sklyanin [4]. In particular, the reflection equation provides the reflection matrices that can be used as the boundary condition for the Bäcklund matrix that one uses to perform the nonlinear mirror image method.

Application of the nonlinear mirror image method to the construction of solutions
In this section, we draw on the results of the previous section about the folding procedure associated to Bäcklund transformations in order to derive explicit solutions of the equations of motions with our integrable boundary conditions. To do so, we first need to review the ISM for the AL model (in our normalisation) on the full line. Then, the construction of solutions on the half-line is implemented as a special Z 2 reduction realised by the Bäcklund matrix B(j, t, z) derived above. For the sake of clarity, we focus on case α = β = 1 2 and γ = −1 in all this section. Then, one gets from now on J = I, τ (z) = 1/z and ω(z) = 1 2 (z − 1/z) 2 . After presenting the main results we need, we will further restrict our attention to the discrete NLS case R j = νQ * j , ν = ±1. As we are interested in soliton solutions, we will set ν = −1.

Review of the ISM for AL on the full lattice in the normalisation (2.3)
In [20], the ISM for the AL system in the normalisation (2.3) was presented, with an emphasis on the Riemann-Hilbert formulation of the inverse part. For our purposes, i.e. the construction of explicit multisoliton solutions, the detailed results of [25] in the original normalisation will be more useful. To be able to use them, we first need to make the connection between the two normalisations as far as ISM is concerned. The best way is to realise that they area related by a gauge transformation as follows.
Let us work with the fields Q j , R j , j ∈ Z. We consider the auxiliary problem on the full lattice (j ∈ Z) in the new normalisation, 4) and, We also consider the auxiliary problem on the full line in the old normalisation, for j ∈ Z, where U (j, z) = (Z + W j ) , (6.8) One can check that the two are related by a gauge transformation of the following form (up to a possible normalisation constant, see below) where F is a scalar function satisfying z)) . (6.11) We fix F by normalising it to 1 as j → −∞ to get We will use this gauge transformation below to transfer the results of [25] to our setup. Let us now review ISM for AL. The initial sequences Q j | t=0 , R j | t=0 are assumed to have finite 1-norm where This ensures the desired analyticity properties below [25]. Also, they are assumed to be such that N j is defined and nonzero for all j ∈ Z. Define and consider the two fundamental solutions Ψ old/new ± (j, t, z) normalised to I as j → ±∞: The gauge transformation between the two Lax pairs implies that where It can be shown that F ∞ is time-independent and hence is just a constant number. As a consequence, we find This is the key relation we need to relate the results of [25] with our setup. One convenient consequence of the new normalisation of L (j, z) is that This is of course consistent with the known fact [25] that det S old (z) = F 2 ∞ for |z| = 1. The analyticity properties in z of the column vectors of Φ old/new ± (j, t, z) and that of the entries of S old/new (z) are crucial for the implementation of the ISM. We see that the gauge transformation (6.17) and relation (6.19) have no consequence on these properties (domain of analyticity, location of zeros). The only consequences are 22 to produce an overall factor F ∞ between the so-called norming constants in the old and new normalisation and to change the normalisation of the scattering coefficients which now reads (see [20]) The change of normalisation of the norming constants can always be absorbed by redefining them. Therefore, in the following, we simply review the results of ISM that we need from [25] and assume that the overall constant F ∞ has been absorbed in the norming constants. In particular, we now focus on the new normalisation and drop the superscript new for conciseness. Let us split the Jost solutions into column vectors and write One shows that Φ L + (j, t, z), Φ R − (j, t, z) and s 22 (z) are analytic functions of z for |z| < 1 and continuous for |z| ≤ 1 while Φ L − (j, t, z), Φ R + (j, t, z) and s 11 (z) are analytic functions of z for |z| > 1 and continuous for |z| ≥ 1. We consider a finite number of simple zeros z n , n = 1, . . . , P of s 22 (z) in the region |z| > 1 and a finite number of zeros z n , n = 1, . . . , P of s 11 (z) in the region |z| < 1. At those zeros, the column vectors of Φ ± (j, t, z) are related by the so-called norming constants as follows It will be convenient to introduce another set of norming constants defined by . (6.26) Equipped with the scattering data, one can implement the inverse part of the method to obtain reconstruction formulas for the solution of the AL system on the full lattice. There is an inherent symmetry on the discrete data. It comes from the fact that the Lax matrices L (j, z) and A (j, z) satisfy the following symmetry relation and As a consequence, the zeros of s 11 (z) (resp. s 22 (z)) come in pairs ±z n and hence P = 2κ for some integer κ ≥ 0 (resp. ±z n , P = 2κ). One can show that the norming constant C n (resp. C n ) associated to z n (resp. z n ) is equal to the norming constant associated to −z n (resp. −z n ). Equations (3.2.102)-(3.2.104) of [25] exploit this symmetry and can be used to derive a nice compact formula for the pure soliton case associated 23 to (half of) the discrete data {z 1 , . . . , z κ ; C 1 , . . . , C κ ; z 1 , . . . , z κ ; C 1 , . . . , C κ }. After some calculations, we find 4 where the κ × κ matrix µ j (t) has the following entries , n, l = 1, . . . , κ , (6.32) and the κ × κ matrix µ j (t) has the following entries , n, l = 1, . . . , κ . (6.33) From now on, we set R j = −Q * j . This leads to an additional symmetry on the Lax pair which in turns implies a symmetry on the scattering data which should be implemented in the above formulas. In short we have, and, on the discrete data, κ = κ and z n = 1 z * n and C n = C * n (z * n ) 2 . (6.37) With these, we can restrict our attention to (6.30) and (6.33). Under this reduction, we can use the trace formulas (3.2.89) of [25] in the pure soliton case to obtain the following useful explicit form of s 11 (z), s 22 (z), suitably normalised to our setup (see in particular (6.21)-(6.22)) (6.38) Since det S(z) = 1 when |z| = 1 and F ∞ > 0 in the DNLS reduction with ν = −1, this fixes the value of F ∞ completely in terms of the scattering data as 5 (6.39)

Symmetry from the mirror image procedure
The important idea behind the mirror image method is that one can obtain solutions of an integrable PDE on the half-line with certain (integrable) boundary conditions by considering solutions of the full line problem associated to scattering data with a special symmetry. One way to obtain the desired symmetry on the scattering data is by interpreting the boundary conditions as arising from a well chosen Bäcklund transformation relating a solution to its "mirror image". This is what we have done above by constructing B(j, t, z) and we use it now to derive the symmetry of the scattering data. We first consider the continuous scattering data Proposition 6.1. The following relation holds on the Jost solutions As a consequence, the scattering matrix satisfies Proof. Let us denote by Φ ± (j, t, z) the Jost solutions of the auxiliary problem with the mirror image solutions. Φ ± (j, t, z) and B(j, t, z) Φ ± (j, t, z) are solutions of the same auxiliary problem. Then, there exist matrices N ± (z) independent of time and position such that It is easy to see that the normalisation in the previous relation is the identity by taking the limit j → ∞. Relations (6.44) and (6.45) lead to relation (6.40) of the proposition. Since Φ − (j, t, z) = Φ + (j, t, z) S(z) and by using (6.40) by replacing j → −j and z → 1/z, we obtain From this and (5.39), one gets (6.41) as desired.
As an important consequence of the Proposition 6.1, one gets We turn to the discrete data. The symmetry relation (6.47) implies that κ = κ and without loss of generality z n = 1 z n , n = 1, . . . , κ . (6.50) Note that the reduction to κ = κ would hold even if we did not consider the DNLS reduction of the AL system. It is the result of the folding symmetry only. The same holds true for the folding symmetry on the scattering data that we discuss below. That being said, since the focus of this section is on DNLS, let us examine the effect of the DNLS reduction on the boundary data before we proceed to showing the effect of the folding symmetry on the discrete data. In particular, f 1 ∞ ∈ R. Proof: The reduction also applies to k − (z) and B(+∞, z) which must therefore satisfy the symmetry condition (6.34). Direct calculation yields the stated results. Proposition 6.2. Under the folding reduction, the discrete data satisfies z n = 1 z n , C n C n = − ϕ(1/z n ) (z n s 11 (z n )) 2 ϕ(z n ) , n = 1, . . . , κ , (6.53) where ϕ(z) is given by (see (5.36), taking into account the DNLS reduction)

Reflected solitons for the DNLS
The three different symmetries on the scattering data (the one inherent to AL, the one coming from the DNLS reduction and the one associated to the folding reduction) are compatible, as we have seen. When they are imposed simultaneously, one can construct solutions of DNLS on the half-lattice with integrable boundary conditions using the solution on the full lattice and restricting to positive integers. More precisely, when all three symmetries hold, the numbers of zeros P of s 22 is P = 2κ where κ is itself an even integer, κ = 2k, as a consequence of the folding symmetry. In that case, we also know that P = P = 4k. Therefore, the discrete data come in octets that produce one soliton each. Each octet is completely determined by one zero z n ∈ C, |z n | > 1, and one norming constant C n ∈ C which are paired as follows, for n = 1, . . . , k, ) , (6.63) ) , (6.64) (1/z n , − ϕ(1/z n ) C n (z n s 11 (z n )) 2 ϕ(z n ) ) , (6.65) ) , (6.66) The first four zeros in each octet correspond to s 11 (z) while the last four correspond to s 22 (z). The explicit formulas (6.38) now take the form The fact that each zero and its opposite have the same norming constant is intrinsic to the model as already discussed, and this has already been taken into account when deriving the explicit formulas (6.30)-(6.33).
The appearance of the discrete data in octets for the AL model with integrable boundary conditions was first shown in [26] and then in [21] in the case of Robin boundary conditions. In our case, the entirety of the effect of our time-dependent boundary conditions is encoded in the function ϕ(z) or, alternatively, in the function f (z) in (6.49). This is in line with the results of [21] for the Robin case. In fact, we can reproduce the known Robin case by choosing c = d = 0 and for instance a = −1, b = 1/2χ: our function f (z) in (6.49) becomes which consistently reproduces the function f (z) of [21] (eq.(7.12)) where our (f 1 ∞ ) 2 plays the role of p ∞ in that paper. In our more general case, a thorough analysis of all the possible values that f 1 ∞ can take among the roots of (6.54), analogously to Corollary 6.4 of [21], would be required to classify all the possible scenarios of solutions that one can construct using the mirror image method. We do not perform this analysis here as the number of cases to consider is much larger than in the Robin case. This technical point does not affect the significance of the results we have obtained. In particular, in the explicit examples to follow, we simply fix compatible numerical values of f 1 ∞ and a, b, c, d to produce plots of solitons being reflected by our boundary conditions.
In Figure 1, we present such a reflected solution in the one-soliton case. Left plots show a one-soliton being reflecting off the boundary at x = −1. We allowed x to be real-valued but highlighted integer values in solid black curves. Right plots show contour plots of the one-soliton being reflected as well as the image soliton on the other side of the boundary (the black vertical line). For comparison, we display 3 types of boundary conditions: our time-dependent case and two particular cases of it which were previously known (Robin and Dirichlet).

Conclusions and outlook
The main message illustrated in this paper is that the Hamiltonian approach and the zero curvature approach to integrable boundary conditions are not separated topics but two faces of the same coin. This interplay is well-known for problems on the line but somehow was not as clear in the case with boundaries, despite the seminal work of Sklyanin [4]. The two aspects inform each other, as we have shown in detail with the AL model in this paper. For instance, the choice of the value of the Bäcklund matrix B(j, t, z) at j = 0 is dictated by the solution of the reflection equation that we consider. Taking non-diagonal solutions of the reflection equation, we showed that new, time-dependent, boundary conditions arise. For the first time, we then developed the nonlinear mirror image method for such boundary conditions and constructed explicit soliton solutions. A study of a continuous model, such as the nonlinear Schrödinger (NLS) equation, along the lines of the present work would be desirable to investigate what kind of more general integrable boundary conditions than the standard Robin boundary conditions can be imposed. This is currently under investigation. We note that the Hamiltonian aspects of this question, for the (vector) NLS equation, have already been discussed in [27]. However, the full connection to zero curvature representation and to the nonlinear mirror method for the construction of solutions remains an open problem.
The results illustrated here rely on the general procedure described in [16] which only involves fundamental features of integrable systems, such as r-matrix and Lax matrix. Therefore, there is no obstacle in principle to apply the same ideas for other models, with the important proviso that the model allows for a natural folding (Z 2 ) symmetry in order to implement the nonlinear mirror image method. The problem of understanding the analog of the mirror image method for models that do not possess a natural Z 2 symmetry is completely open and rather fascinating. So far, the only alternative to discuss such models (e.g. KdV) on the half-line with boundary conditions is to use the so-called unified transform [28]. It is a completely open problem to investigate integrable time-dependent boundary conditions of the kind we found in that setup.