# Bäcklund transformations for many-body systems related to KdV

## 1. Introduction

Bäcklund transformations (BTs) are an important aspect of the theory of integrable systems which have traditionally been studied in the context of evolution equations. However, more recently there has been much interest in discrete systems or integrable mappings . Within the modern approach to separation of variables (reviewed by Sklyanin in ) this has led to the study of BTs for finite-dimensional Hamiltonian systems . The latter are canonical transformations including a Bäcklund parameter $`\lambda `$, and apart from being interesting integrable mappings in their own right they also lead to separation of variables when $`n`$ such mappings are applied to an integrable system with $`n`$ degrees of freedom. The sequence of Bäcklund parameters $`\lambda _j`$ together with a set of conjugate variables $`\mu _j`$ constitute the separation variables, and satisfy a new property called spectrality introduced in .

We proceed to develop these ideas with some new examples of BTs for $`n`$-body systems, namely the many-body generalisation of the case (ii) integrable Hénon-Heiles system, the Garnier system and the Neumann system on the sphere (see ). It is known that the case (ii) Hénon-Heiles system is equivalent to the stationary flow of fifth-order KdV , while the Garnier and Neumann systems may be obtained as restricted flows of the KdV hierarchy . Thus we derive BTs for these systems by reduction of the standard BT for KdV, which arises from the Darboux-Crum transformation  for Schrödinger operators. The restriction of the Darboux transformation to the stationary flows of the modified (mKdV) hierarchy has been discussed in .

In the following section we describe how the reduction works in general, before specialising these considerations to each particular system and presenting the associated generating function for the BT. We note that these systems are examples of the reduced Gaudin magnet , so that we have the following Lax matrix

(1.1) 
$$L(u)=\underset{j=1}{\overset{n}{}}\frac{\mathrm{}_j}{ua_j}+B(u),\mathrm{}_j=\left(\begin{array}{cc}S_j^3& S_j^{}\\ S_j^+& S_j^3\end{array}\right)$$

where (up to scaling) the $`S_j`$ satisfy $`n`$ independent copies of the standard $`sl(2)`$ algebra:

(1.2) 
$$\{S_j^3,S_k^\pm \}=\pm 2\delta _{jk}S_k^\pm ,\{S_j^+,S_k^{}\}=4\delta _{jk}S_k^3.$$

For the Hénon-Heiles and Garnier systems the matrix $`B(u)`$ is respectively quadratic and linear in the spectral parameter $`u`$, while for the Neumann system it is independent of $`u`$ and turns out to be constant due to the constraint that the particles lie on the sphere (hence the Poisson algebra (1.2) must be modified by Dirac reduction).

We have constructed the BT for the (non-reduced) $`sl(2)`$ Gaudin magnet with quasi-periodic boundary condition in , while some preliminary results on the classical Garnier system and two-body Hénon-Heiles system first appeared in .

## 2. Classical integrable systems and KdV

### 2.1. Restricting the BT

As is well known, the Darboux-Crum transformation  consists of mapping the Schrödinger operator $`_t^2+V\lambda `$ to another operator $`_t^2+\stackrel{~}{V}\lambda `$ by factorizing the former and then reversing the order of factorisation. Given an eigenfunction $`\varphi `$ satisfying

$$(_t^2+V\lambda )\varphi =0$$

we may set $`y=(\mathrm{log}[\varphi ])_t`$ and then

(2.1) 
$$V=y_ty^2+\lambda ,\stackrel{~}{V}=y_ty^2+\lambda ;$$

for $`\lambda =0`$ this is just the Miura map for KdV. Also given another eigenfunction $`\psi `$ of the Schrödinger operator with potential $`V`$ for a different spectral parameter $`u`$ we have

$$(_t^2+Vu)\psi =0,(_t^2+\stackrel{~}{V}u)\stackrel{~}{\psi }=0$$

where the transformation to the new eigenfunction $`\stackrel{~}{\psi }`$ and its derivative may be given in matrix form as

(2.2) 
$$\left(\begin{array}{c}\stackrel{~}{\psi }_t\\ \stackrel{~}{\psi }\end{array}\right)=k\left(\begin{array}{cc}y& y^2+u\lambda \\ 1& y\end{array}\right)\left(\begin{array}{c}\psi _t\\ \psi \end{array}\right)$$

for any constant $`k`$. From (2.1) follows the standard formula for the Darboux-Bäcklund transformation of KdV, $`\stackrel{~}{V}=V+2(\mathrm{log}[\varphi ])_{tt}`$.

For what follows it will also be necessary to consider a product of eigenfunctions for a Schrödinger operator with potential $`V`$ and eigenvalue $`u`$,

$$f=\psi \psi ^{}$$

with Wronskian $`\psi _t\psi ^{}\psi \psi _t^{}=2m`$. It is well known that $`f`$ satisfies the Ermakov-Pinney equation 

(2.3) 
$$ff_{tt}\frac{1}{2}f_t^2+2(Vu)f^2+2m^2=0.$$

If we now transform $`\psi `$ and $`\psi ^{}`$ according to (2.2) then we find a new product of eigenfunctions $`\stackrel{~}{f}=\stackrel{~}{\psi }\stackrel{~}{\psi }^{}`$ satisfying the same Ermakov-Pinney equation but with $`V`$ replaced by $`\stackrel{~}{V}`$, given explicitly by

(2.4) 
$$\stackrel{~}{f}=(\lambda u)^1\frac{(Z^2m^2)}{f},Z=\frac{1}{2}f_tyf,$$

where we have set $`k^2=(\lambda u)^1`$ to ensure that the transformed eigenfunctions have the same Wronskian $`2m`$. It is also straightforward to show that, in terms of $`\stackrel{~}{f}`$, the quantity $`Z`$ can be written as $`Z=\frac{1}{2}\stackrel{~}{f}_ty\stackrel{~}{f}`$ (see ).

We can now describe how this transformation restricts to the finite-dimensional Hamiltonian systems presented below. The systems are expressed in variables $`(q_j,p_j)`$ which appear in the Lax matrix (1.1) via the identification 

$$S_j^3=p_jq_j,S_j^{}=p_j^2+\frac{m_j^2}{q_j^2},S_j^+=q_j^2.$$

For Hénon-Heiles and Garnier the non-vanishing Poisson brackets are the standard ones $`\{p_j,q_k\}=\delta _{jk}`$ which provide a realization of the algebra (1.2); for the Neumann system on the sphere the bracket must be modified by Dirac reduction.

All of the systems are Liouville integrable, and thus have a complete set of Hamiltonians in involution, but for these purposes we concentrate on the Hamiltonian $`h`$ generating the flow corresponding to $`t`$ above (in KdV theory this is usually denoted $`x`$, the spatial variable). For this flow the Lax equation $`L_t=[N,L]`$ is the compatibility condition for the linear system

(2.5) 
$$L(u)\mathrm{\Psi }=v\mathrm{\Psi },\mathrm{\Psi }_t=N\mathrm{\Psi };N=\left(\begin{array}{cc}0& uV(q_j,p_j)\\ 1& 0\end{array}\right).$$

Observe that the second part of the linear system is just a Schrödinger equation for the potential $`V`$; for Neumann and Garnier this is a function of $`(q_j,p_j)`$ for $`j=1,\mathrm{},n`$, while for Hénon-Heiles there is an extra pair of conjugate variables $`(q_{n+1},p_{n+1})`$ such that $`Vq_{n+1}`$.

The equations of motion generated by this Hamiltonian take the form $`q_{j,t}=p_j`$ and

(2.6) 
$$p_{j,t}=q_{j,tt}=(a_jV(q_k,p_k))q_j\frac{m_j^2}{q_j^3}$$

for $`j=1,\mathrm{},n`$; for Hénon-Heiles there are also equations for $`q_{n+1}`$ and $`p_{n+1}=q_{n+1,t}`$. The important thing to observe is that (2.6) is equivalent to the fact that $`S_j^+=q_j^2`$ satisfies the Ermakov-Pinney equation (2.3) corresponding to a Schrödinger equation with potential $`V`$ and eigenvalue $`a_j`$. Thus to obtain a Bäcklund transformation for these many-body systems we simply apply a Darboux-Crum transformation (2.1) to the potential $`V=V(q_j,p_j)`$ to obtain $`\stackrel{~}{V}=V(\stackrel{~}{q}_j,\stackrel{~}{p}_j)`$, and then we know that the solutions of the Ermakov-Pinney equation must transform according to (2.4). By this procedure we may explicitly construct the BT for the many-body systems below (or for any restricted flow of KdV), and it is then simple to calculate the generating function $`F(q_j,\stackrel{~}{q}_j)`$ of this canonical transformation, such that

$$dF=\underset{j}{}(p_jdq_j\stackrel{~}{p}_jd\stackrel{~}{q}_j).$$

The discrete Lax equation for the BT,

$$ML=\stackrel{~}{L}M$$

where $`\stackrel{~}{L}=L(\stackrel{~}{q_j},\stackrel{~}{p_j};u)`$, is necessary to ensure the preservation of the spectral curve $`det(vL(u))=0`$ (so that all the Hamiltonians in involution are preserved). This follows immediately from the properties of the Darboux-Crum transformation, since we know that the vector $`\mathrm{\Psi }`$ in the linear system (2.5) must transform according to (2.2), and hence we may take (setting $`k=1`$)

(2.7) 
$$M=\left(\begin{array}{cc}y& y^2+u\lambda \\ 1& y\end{array}\right).$$

Of course we must determine $`y`$ as a function of the dynamical variables. In the Garnier and Hénon-Heiles cases it turns out that the potential depends on coordinates only, $`V=V(q_j)`$, and so by adding the two equations in (2.1) we obtain

$$y(q_j,\stackrel{~}{q_j})=\pm \sqrt{\lambda \frac{1}{2}(V+\stackrel{~}{V})};$$

to obtain the correct continuum limit of the discrete dynamics it is necessary to take the negative branch of the square root (see ). For the Neumann system $`V`$ depends on both coordinates and momenta, so the above does not yield $`y(q_j,\stackrel{~}{q_j})`$.

There is another way of writing $`L`$ which arises more naturally via reduction from the zero curvature representation of the KdV hierarchy , viz

$$L(u)=\left(\begin{array}{cc}\frac{1}{2}\mathrm{\Pi }_t& \hfill \frac{1}{2}\mathrm{\Pi }_{tt}+(uV)\mathrm{\Pi }\\ \mathrm{\Pi }& \hfill \frac{1}{2}\mathrm{\Pi }_t\end{array}\right)$$

where

(2.8) 
$$\mathrm{\Pi }(u)=\underset{j=1}{\overset{n}{}}\frac{q_j^2}{ua_j}+\mathrm{\Delta }(u).$$

$`\mathrm{\Delta }`$ is a polynomial in $`u`$ fixing the dynamical term $`B`$ in (1.1); we shall present the appropriate $`\mathrm{\Delta }`$ and $`B`$ in each case below. Clearly the $`t`$ derivatives of $`\mathrm{\Pi }`$ can be rewritten using the equations of motion to yield (1.1).

Finally if we write the (hyper-elliptic) spectral curve as

$$v^2=R(u)$$

then it is easy to check that the spectrality property  is satisfied for these systems, in the sense that defining the conjugate variable to $`\lambda `$ by

$$\mu =2\frac{F}{\lambda }$$

we find that

$$L(\lambda )\mathrm{\Omega }=\mu \mathrm{\Omega }$$

with eigenvector $`\mathrm{\Omega }=(y,1)^T`$, or in other words $`\mu ^2=R(\lambda )`$ so that $`(\lambda ,\mu )`$ is a point on the spectral curve. Note that (as for the examples in ) this eigenvector spans the kernel of $`M`$,

$$M(\lambda )\mathrm{\Omega }=0.$$

We can also write $`y`$ explicitly in terms of both the old and the new variables related by the BT, thus:

(2.9) 
$$y(q_j,p_j)=\frac{\mathrm{\Pi }_t(\lambda )+2\mu }{2\mathrm{\Pi }(\lambda )},y(\stackrel{~}{q_j},\stackrel{~}{p_j})=\frac{(\stackrel{~}{\mathrm{\Pi }}_t(\lambda )2\mu )}{2\stackrel{~}{\mathrm{\Pi }}(\lambda )};$$

clearly we denote $`\stackrel{~}{\mathrm{\Pi }}(\lambda )=\mathrm{\Pi }(\stackrel{~}{q_j},\stackrel{~}{p_j};\lambda )`$.

### 2.2. Generalised Hénon-Heiles system

For the many-body generalisation of case (ii) integrable Hénon-Heiles system, the Hamiltonian generating the $`t`$ flow takes the form

$$h=\frac{1}{2}\underset{j=1}{\overset{n+1}{}}p_j^2+q_{n+1}^3+q_{n+1}\left(\frac{1}{2}\underset{j=1}{\overset{n}{}}q_j^2+c\right)\frac{1}{2}\underset{j=1}{\overset{n}{}}\left(a_jq_j^2+\frac{m_j^2}{q_j^2}\right).$$

The original case (ii) integrable Hénon-Heiles system corresponds to $`n=1`$ with $`c=m_j=a_j=0`$. The link between stationary fifth-order KdV and the type (ii) system was noted by Fordy in , although this was anticipated in work of Weiss , who used Painlevé analysis to derive a BT and associated linear problem (a similar result also appears in ). None of these authors wrote a BT explicitly as a canonical transformation with parameter, although (without parameter) this was done for a non-autonomous version in .

For the Lax matrix $`L`$ of the generalised $`(n+1)`$-body Hénon-Heiles system we have $`\mathrm{\Delta }=16u8q_{n+1}`$ so that the extra term $`B(u)`$ is given by

$$B=\left(\begin{array}{cc}4p_{n+1}& E\\ 16u8q_{n+1}& 4p_{n+1}\end{array}\right),E=16u^2+8q_{n+1}u4q_{n+1}^2\underset{j=1}{\overset{n}{}}q_j^24c.$$

The equations of motion for $`h`$ imply that the squares of the first $`n`$ coordinates $`q_j^2`$ satisfy the Ermakov-Pinney equation (2.3) for $`m=m_j`$ with

$$V=q_{n+1}$$

and eigenvalue $`a_j`$. Thus the BT for the system can be calculated directly by applying the Darboux-Crum transformation to $`V=q_{n+1}`$, to yield $`\stackrel{~}{V}=\stackrel{~}{q}_{n+1}`$, and applying (2.4) to each $`q_j^2`$ for $`j=1,\mathrm{},n`$.

After some calculation the generating function for this canonical transformation is found to be

$$F(q_j,\stackrel{~}{q}_j;\lambda )=\underset{j=1}{\overset{n}{}}\left(Z_j+\frac{m_j}{2}\mathrm{log}\left[\frac{Z_jm_j}{Z_j+m_j}\right]\right)+\frac{16}{5}y^5+4(q_{n+1}+\stackrel{~}{q}_{n+1})y^3$$

$$+\left(2q_{n+1}^2+2q_{n+1}\stackrel{~}{q}_{n+1}+2\stackrel{~}{q}_{n+1}^2+\frac{1}{2}\underset{j=1}{\overset{n}{}}(q_j^2+\stackrel{~}{q}_j^2)+2c\right)y,$$

where we have found it convenient to use the quantities $`Z_j(q_j,\stackrel{~}{q}_j)`$ and $`y(q_j,\stackrel{~}{q}_j)`$ defined by

(2.10) 
$$Z_j^2=m_j^2+(\lambda a_j)q_j^2\stackrel{~}{q}_j^2,$$

and

$$y=\sqrt{\lambda \frac{1}{2}(q_{n+1}+\stackrel{~}{q}_{n+1})}.$$

In order to check the spectrality property, we have explicitly found that the eigenvalue of $`L(\lambda )`$ with eigenvector $`\mathrm{\Omega }=(y,1)^T`$ can be written as

$$\mu (q_j,\stackrel{~}{q}_j;\lambda )=\underset{j=1}{\overset{n}{}}\frac{Z_j}{\lambda a_j}\frac{1}{y}\frac{F}{y},$$

which precisely equals $`2\frac{F}{\lambda }`$ as required.

### 2.3. Garnier system

For the Garnier system the $`t`$ flow is generated by the Hamiltonian

$$h=\frac{1}{2}\underset{j=1}{\overset{n}{}}p_j^2+\frac{1}{2}(\underset{j=1}{\overset{n}{}}q_j^2)^2\frac{1}{2}\underset{j=1}{\overset{n}{}}\left(a_jq_j^2+\frac{m_j^2}{q_j^2}\right).$$

This differs from the traditional Garnier system as in  by the inclusion of extra inverse square terms. The Newton equations for the $`q_j`$ are

$$q_{j,tt}+2(\underset{k}{}q_k^2)q_j=a_jq_j\frac{m_j^2}{q_j^3},$$

so clearly for the standard restricted flows of KdV , when $`m_j=0`$, each $`q_j`$ is an eigenfunction of a Schrödinger operator with potential

$$V=2\underset{j}{}q_j^2$$

and eigenvalue $`a_j`$, while in general $`q_j^2`$ is a product of eigenfunctions satisfying the Ermakov-Pinney equation  for $`m=m_j`$.

The Lax matrix of the Garnier system has $`\mathrm{\Delta }=1`$, so $`L`$ takes the form (1.1) with

$$B=\left(\begin{array}{cc}0& u_jq_j^2\\ 1& 0\end{array}\right).$$

Applying the Darboux-Crum transformation we obtain a new potential

$$\stackrel{~}{V}=2\underset{j}{}\stackrel{~}{q}_j^2,$$

and the corresponding BT induced on the Garnier system is equivalent to gauging $`L`$ by the matrix $`M`$ of the form (2.7) with

$$y=\sqrt{\lambda \underset{j}{}(q_j^2+\stackrel{~}{q}_j^2)}.$$

Finally we can calculate the generating function for this BT, which may be written as follows:

$$F(q_j,\stackrel{~}{q}_j;\lambda )=\underset{j=1}{\overset{n}{}}\left(Z_j+\frac{m_j}{2}\mathrm{log}\left[\frac{Z_jm_j}{Z_j+m_j}\right]\right)\frac{1}{3}y^3,$$

where $`y(q_j,\stackrel{~}{q}_j)`$ is as above and $`Z_j`$ is given by the same expression (2.10) as for Hénon-Heiles. In  we derived this generating function for the special case $`m_j=0`$ when the logarithm terms do not appear. To check spectrality we notice that $`L(\lambda )`$ has eigenvalue

$$\mu (q_j,\stackrel{~}{q}_j;\lambda )=\underset{j=1}{\overset{n}{}}\frac{Z_j}{\lambda a_j}+y$$

with eigenvector $`\mathrm{\Omega }`$, and so we see that $`\mu =2\frac{F}{\lambda }`$.

### 2.4. Neumann system on the sphere

For the Neumann system the $`t`$ flow is generated by

$$h=\frac{1}{2}\underset{j=1}{\overset{n}{}}p_j^2\frac{1}{2}\underset{j=1}{\overset{n}{}}\left(a_jq_j^2+\frac{m_j^2}{q_j^2}\right).$$

Once again this has extra inverse square terms compared with the standard Neumann system . The Poisson bracket for this system is modified by constraining the particles to lie on a sphere, so that

(2.11) 
$$(q,q)\underset{j}{}q_j^2=const,(q,p)\underset{j}{}q_jp_j=0$$

which results in the non-vanishing Dirac brackets

(2.12) 
$$\{p_j,q_k\}=\delta _{jk}\frac{q_jq_k}{(q,q)},\{p_j,p_k\}=\frac{q_jp_kq_kp_j}{(q,q)}.$$

With this bracket the Hamilton equations are $`q_{j,t}=p_j`$ and (2.6) with

$$V=(q,q)^1\underset{j}{}\left(p_j^2+a_jq_j^2\frac{m_j^2}{q_j^2}\right).$$

The Lax matrix for the Neumann system arises by setting $`\mathrm{\Delta }=0`$, which in (1.1) gives the following matrix $`B`$:

$$B=\left(\begin{array}{cc}0& (q,q)\\ 0& 0\end{array}\right).$$

In fact if we start from the linear system (2.5) and leave $`V`$ unspecified then (2.6) as well as the constraint $`(q,q)_t=0`$ are consequences of the Lax equation, and together these are sufficient to determine the form of $`V`$; this is also how the equations for the constrained Neumann system arise in a Lagrangian approach .

Given that the phase space is now degenerate with two Casimirs given by (2.11), it would appear that the standard sort of generating function will no longer be appropriate for describing a BT. It turns out that we can apply the Darboux-Crum transformation as before, and transform the quantities $`q_j^2`$ according to (2.4). In this way we obtain new variables $`\stackrel{~}{q}_j(q_k,p_k)`$ and $`\stackrel{~}{p}_j(q_k,p_k)`$, which are naturally written with the use of the quantity $`y(q_k,p_k)`$ given by the first formula in (2.9); on the Lax matrix this transformation arises by gauging with $`M`$ as in (2.7). Similarly the transformation can be inverted to give $`q_j(\stackrel{~}{q}_k,\stackrel{~}{p}_k)`$ and $`p_j(\stackrel{~}{q}_k,\stackrel{~}{p}_k)`$ written in terms of $`y(\stackrel{~}{q}_k,\stackrel{~}{p}_k)`$ given by the right hand formula of (2.9).

However, it would still be nice to write a generating function for this transformation. We have found that if we formally take

$$F(q_j,\stackrel{~}{q}_j;\lambda )=\underset{j=1}{\overset{n}{}}\left(Z_j+\frac{m_j}{2}\mathrm{log}\left[\frac{Z_jm_j}{Z_j+m_j}\right]+\frac{1}{2}y(q_j^2\stackrel{~}{q}_j^2)\right)$$

with $`Z_j`$ given by (2.10) as usual, and regard $`y`$ as a sort of Lagrange multiplier (independent of the coordinates and $`\lambda `$), then we do indeed obtain the correct expressions

$$p_j=\frac{F}{q_j},\stackrel{~}{p}_j=\frac{F}{\stackrel{~}{q}_j},$$

but these contain $`y`$ which is unspecified. If we then require that the constraints (2.11) are preserved under the BT applied from old to new variables or vice-versa, then in either direction the constraints are preserved if and only if $`y`$ satisfies a quadratic equation with solution given respectively by the formulae (2.9). Alternatively if we require spectrality then second component of the equation $`L(\lambda )\mathrm{\Omega }=\mu \mathrm{\Omega }`$ gives

$$\mu (q_j,\stackrel{~}{q}_j;\lambda )=\underset{j=1}{\overset{n}{}}\frac{Z_j}{\lambda a_j}=2\frac{F}{\lambda }$$

as required, while the first component gives (after making use of the formula (2.10) and the BT)

$$\mu =\underset{j=1}{\overset{n}{}}\frac{Z_j}{\lambda a_j}+\frac{1}{y}\underset{j}{}(q_j^2\stackrel{~}{q}_j^2).$$

Hence spectrality requires that the second term vanishes, and so the first constraint (2.11) is preserved; the preservation of the second constraint is then an algebraic consequence of the BT.

Thus we see that for this BT we can write the new variables as functions of the old and vice-versa, but a formula for $`y(q_j,\stackrel{~}{q}_j;\lambda )`$ is lacking. Also this discretization of the Neumann system is apparently new, since it is exact (preserving the Lax matrix for the continuous system) unlike the Veselov or Ragnisco discretizations discussed in .

## 3. Conclusions

It would also be interesting to look at BTs with parameter in the non-autonomous case , where deformation with respect to the Bäcklund parameter would probably have to be introduced (corresponding to the associated isomonodromy problem).

## 4. Acknowledgements

ANWH thanks the Leverhulme Trust for providing a Study Abroad Studentship in Rome, and is grateful to J. Harnad and Y. Suris for useful conversations. VBK acknowledges the support from the EPSRC and the support from Istituto Nazionale di Fisica Nucleare for his visit to Rome. The authors would also like to thank the organisers of the meeting Integrable Systems: Solutions and Transformations in Guardamar, Spain (June 1998) where some of this work was carried out.