Direct linearisation approach to discrete integrable systems associated with $\mathbb{Z}_\mathcal{N}$ graded Lax pairs

TL;DR

This paper links Fordy-Xenitidis discrete integrable systems to linear integral equations, revealing their solution structure and connecting them to the discrete Gel'fand-Dikii hierarchy, while also presenting their bilinear form and tau functions.

Contribution

It establishes a novel connection between Fordy-Xenitidis systems and linear integral equations, providing explicit solutions and bilinear forms.

Findings

Solution structure of Fordy-Xenitidis systems revealed

Connection to discrete Gel'fand-Dikii hierarchy established

Bilinear form and tau functions derived

Abstract

Fordy and Xenitidis [J. Phys. A: Math. Theor. 50 (2017) 165205] recently proposed a large class of discrete integrable systems which include a number of novel integrable difference equations, from the perspective of $\mathbb{Z}_\mathcal{N}$ graded Lax pairs, without providing solutions. In this paper, we establish the link between the Fordy-Xenitidis discrete systems in coprime case and linear integral equations in certain form, which reveals solution structure of these equations. The bilinear form of the Fordy-Xenitidis integrable difference equations is also presented together with the associated general tau function. Furthermore, the solution structure explains the connections between the Fordy-Xenitidis novel models and the discrete Gel'fand-Dikii hierarchy.

Equations

$$(p+q+u_{n,m}-u_{n+1,m+1})(p-q+u_{n,m+1}-u_{n+1,m})=p^2-q^2.$$

$$p(v_{n,m}{v}_{n,m+1}-v_{n+1,m}{v}_{n+1,m+1})=q(v_{n,m}{v}_{n+1,m}-v_{n,m+1}{v}_{n+1,m+1})$$

$$\frac{(z_{n,m}-z_{n,m+1})(z_{n+1,m}-z_{n+1,m+1})}{(z_{n,m}-z_{n+1,m})(z_{n,m+1}-z_{n+1,m+1})}=\frac{p^2}{q^2}$$

$$(p-q{)}^2{\tau }_{n,m}{\tau }_{n+2,m+2}-(p+q{)}^2{\tau }_{n+2,m}{\tau }_{n,m+2}+4pq{\tau }_{n+1,m+1}^2=0.$$

$$(p+q){\tau }_{n,m+1}{\tau }_{n+2,m}+(p-q){\tau }_{n,m}{\tau }_{n+2,m+1}=2p{\tau }_{n+1,m}{\tau }_{n+1,m+1},$$

$$(p+q){\tau }_{n+1,m}{\tau }_{n,m+2}-(p-q){\tau }_{n,m}{\tau }_{n+1,m+2}=2q{\tau }_{n,m+1}{\tau }_{n+1,m+1}.$$

$$\frac{p^3-q^3}{u_{n+2,m}-u_{n+1,m+1}}-\frac{p^3-q^3}{u_{n+1,m+1}-u_{n,m+2}}$$

$$=u_{n+2,m+2}(u_{n+2,m+1}-u_{n+1,m+2})-u_{n+1,m}{u}_{n+2,m+1}+u_{n,m+1}{u}_{n+1,m+2}+u_{n,m}(u_{n+1,m}-u_{n,m+1}),$$

$$\frac{v_{n,m}}{v_{n+1,m}}-\frac{v_{n,m}}{v_{n,m+1}}+\frac{v_{n+2,m+1}}{v_{n+2,m+2}}-\frac{v_{n+1,m+2}}{v_{n+2,m+2}}$$

$$=\left(\frac{p^3{v}_{n+1,m+1}-q^3{v}_{n,m+2}}{v_{n,m+2}-v_{n+1,m+1}}\right)\frac{v_{n+1,m+2}}{v_{n,m+1}}-\left(\frac{p^3{v}_{n+2,m}-q^3{v}_{n+1,m+1}}{v_{n+1,m+1}-v_{n+2,m}}\right)\frac{v_{n+2,m+1}}{v_{n+1,m}},$$

$$\frac{(z_{n+2,m+2}-z_{n+1,m+2})(z_{n,m+2}-z_{n+1,m+1})(z_{n,m+1}-z_{n,m})}{(z_{n+2,m+2}-z_{n+2,m+1})(z_{n+2,m}-z_{n+1,m+1})(z_{n+1,m}-z_{n,m})}$$

$$=\frac{p^3(z_{n+1,m+2}-z_{n,m+2})(z_{n+1,m+1}-z_{n,m+1})-q^3(z_{n+1,m+2}-z_{n+1,m+1})(z_{n,m+2}-z_{n,m+1})}{q^3(z_{n+2,m+1}-z_{n+2,m})(z_{n+1,m+1}-z_{n+1,m})-p^3(z_{n+2,m+1}-z_{n+1,m+1})(z_{n+2,m}-z_{n+1,m})}.$$

$$\left(p^2+pq+q^2\right)({\tau }_{n+1,m}{\tau }_{n,m+2}{\tau }_{n+2,m+1}+{\tau }_{n,m+1}{\tau }_{n+2,m}{\tau }_{n+1,m+2})$$

$$=3p^2{\tau }_{n+1,m}{\tau }_{n+1,m+1}{\tau }_{n+1,m+2}+3q^2{\tau }_{n,m+1}{\tau }_{n+1,m+1}{\tau }_{n+2,m+1}-(p-q{)}^2{\tau }_{n,m}{\tau }_{n+1,m+1}{\tau }_{n+2,m+2}.$$

$$\mathbf{\Sigma }={\left(\begin{array}{ccccc}0& 1& & & \\ & 0& 1& & \\ & & \ddots & \ddots & \\ & & & \ddots & 1\\ 1& & & & 0\end{array}\right)}_{\mathcal{N}\times \mathcal{N}},$$

$$\frac{p+u_{n,m+1}^{(r+1+\alpha )}-u_{n+1,m+1}^{(r)}}{p+u_{n,m}^{(r+1+\alpha +\beta )}-u_{n+1,m}^{(r+\beta )}}=\frac{q+u_{n+1,m}^{(r+1+\beta )}-u_{n+1,m+1}^{(r)}}{q+u_{n,m}^{(r+1+\alpha +\beta )}-u_{n,m+1}^{(r+\alpha )}},u_{n,m}^{(r+\mathcal{N})}=u_{n,m}^{(r)},$$

$$p\left(\frac{v_{n+1,m+1}^{(r)}}{v_{n,m+1}^{(r+\alpha )}}-\frac{v_{n+1,m}^{(r+1+\beta )}}{v_{n,m}^{(r+1+\alpha +\beta )}}\right)=q\left(\frac{v_{n+1,m+1}^{(r)}}{v_{n+1,m}^{(r+\beta )}}-\frac{v_{n,m+1}^{(r+1+\alpha )}}{v_{n,m}^{(r+1+\alpha +\beta )}}\right),v_{n,m}^{(r+\mathcal{N})}=v_{n,m}^{(r)},$$

$$\prod _{r=0}^{\mathcal{N}-1}\left(p+u_{n,m}^{(r+1+\alpha )}-u_{n+1,m}^{(r)}\right)=p^{\mathcal{N}},\prod _{r=0}^{\mathcal{N}-1}\left(q+u_{n,m}^{(r+1+\beta )}-u_{n,m+1}^{(r)}\right)=q^{\mathcal{N}}and\prod _{r=0}^{\mathcal{N}-1}{v}_{n,m}^{(r)}=1,$$

$$\frac{p+u_{n,m+1}^{(r+1)}-u_{n+1,m+1}^{(r)}}{p+u_{n,m}^{(r+1)}-u_{n+1,m}^{(r)}}=\frac{q+u_{n+1,m}^{(r+1)}-u_{n+1,m+1}^{(r)}}{q+u_{n,m}^{(r+1)}-u_{n,m+1}^{(r)}},r=0,1,\cdots ,\mathcal{N}-3,$$

$$\left(p+q+u_{n,m}^{(0)}-u_{n+1,m+1}^{(\mathcal{N}-2)}\right)\left(p-q+u_{n,m+1}^{(\mathcal{N}-2)}-u_{n+1,m}^{(\mathcal{N}-2)}\right)$$

$$=\frac{p^{\mathcal{N}}}{{\prod }_{r=0}^{\mathcal{N}-3}\left(p+u_{n,m}^{(r+1)}-u_{n+1,m}^{(r)}\right)}-\frac{q^{\mathcal{N}}}{{\prod }_{r=0}^{\mathcal{N}-3}\left(q+u_{n,m}^{(r+1)}-u_{n,m+1}^{(r)}\right)}$$

$$p\left(\frac{v_{n+1,m+1}^{(r)}}{v_{n,m+1}^{(r)}}-\frac{v_{n+1,m}^{(r+1)}}{v_{n,m}^{(r+1)}}\right)=q\left(\frac{v_{n+1,m+1}^{(r)}}{v_{n+1,m}^{(r)}}-\frac{v_{n,m+1}^{(r+1)}}{v_{n,m}^{(r+1)}}\right),r=0,1,\cdots ,\mathcal{N}-3,$$

$$p\left(\frac{v_{n+1,m+1}^{(\mathcal{N}-2)}}{v_{n,m+1}^{(\mathcal{N}-2)}}-\frac{{\prod }_{r=0}^{\mathcal{N}-2}{v}_{n,m}^{(r)}}{{\prod }_{r=0}^{\mathcal{N}-2}{v}_{n+1,m}^{(r)}}\right)=q\left(\frac{v_{n+1,m+1}^{(\mathcal{N}-2)}}{v_{n+1,m}^{(\mathcal{N}-2)}}-\frac{{\prod }_{r=0}^{\mathcal{N}-2}{v}_{n,m}^{(r)}}{{\prod }_{r=0}^{\mathcal{N}-2}{v}_{n,m+1}^{(0)}}\right)$$

$$\frac{p+u_{n,m+1}^{(0)}-u_{n+1,m+1}^{(0)}}{q+u_{n+1,m}^{(0)}-u_{n+1,m+1}^{(0)}}=\frac{p^{\mathcal{N}}{\prod }_{r=0}^{\mathcal{N}-2}\left(q+u_{n,m}^{(r)}-u_{n,m+1}^{(r)}\right)}{q^{\mathcal{N}}{\prod }_{r=0}^{\mathcal{N}-2}\left(p+u_{n,m}^{(r)}-u_{n+1,m}^{(r)}\right)},$$

$$\frac{p+u_{n,m+1}^{(r)}-u_{n+1,m+1}^{(r)}}{q+u_{n+1,m}^{(r)}-u_{n+1,m+1}^{(r)}}=\frac{p+u_{n,m}^{(r-1)}-u_{n+1,m}^{(r-1)}}{q+u_{n,m}^{(r-1)}-u_{n,m+1}^{(r-1)}},r=1,2,\cdots ,\mathcal{N}-3,$$

$$\left(\frac{p{v}_{n+1,m+1}-q{v}_{n+2,m}}{p{v}_{n,m+2}-q{v}_{n+1,m+1}}\right)\left(\frac{p^2{v}_{n+1,m+1}-q^2{v}_{n,m+2}}{p^2{v}_{n+2,m}-q^2{v}_{n+1,m+1}}\right)=\frac{X_{n,m}{X}_{n,m+1}{X}_{n+1,m+1}}{Y_{n,m}{Y}_{n+1,m}{Y}_{n+1,m+1}}$$

$$X_{n,m}\doteq p{v}_{n,m}{v}_{n+1,m}{v}_{n+1,m+1}+qand{Y}_{n,m}\doteq q{v}_{n,m}{v}_{n,m+1}{v}_{n+1,m+1}+p,$$

$${\varphi }_{n+1,m}^{(r)}=\left(p+u_{n,m}^{(r+1+\alpha )}-u_{n+1,m}^{(r)}\right){\varphi }_{n,m}^{(r+\alpha )}+{\varphi }_{n,m}^{(r+1+\alpha )},$$

$${\varphi }_{n,m+1}^{(r)}=\left(q+u_{n,m}^{(r+1+\beta )}-u_{n,m+1}^{(r)}\right){\varphi }_{n,m}^{(r+\beta )}+{\varphi }_{n,m}^{(r+1+\beta )},$$

$$p+u_{n,m}^{(r+1+\alpha )}-u_{n+1,m}^{(r)}=p\frac{v_{n+1,m}^{(r)}}{v_{n,m}^{(r+\alpha )}}andq+u_{n,m}^{(r+1+\beta )}-u_{n,m+1}^{(r)}=q\frac{v_{n,m+1}^{(r)}}{v_{n,m}^{(r+\beta )}}.$$

$${\varphi }_{n,m}^{(r)}(k)+\sum _{j\in J}{\int }_{{\Gamma }_j}\mathrm{d}{\lambda }_j(l)\frac{{\rho }_{n,m}^{(r)}(k){\sigma }_{n,m}^{(r)}({\omega }^{j}l)}{k-{\omega }^{j}l}{\varphi }_{n,m}^{(r)}(l)={\rho }_{n,m}^{(r)}(k),$$

$${\rho }_{n,m}^{(r)}(k)=(k^{\alpha }(p+k){)}^n(k^{\beta }(q+k){)}^m{k}^{r}and{\sigma }_{n,m}^{(r)}(k^{\prime })=(k^{\prime \alpha }(p+k^{\prime }){)}^{-n}(k^{\prime \beta }(q+k^{\prime }){)}^{-m}{k}^{\prime -r};$$

$$u_{n,m}^{(r)}\doteq \sum _{j\in J}{\int }_{{\Gamma }_j}\mathrm{d}{\lambda }_j(k){\varphi }_{n,m}^{(r)}(k){\sigma }_{n,m}^{(r)}({\omega }^{j}k)and{v}_{n,m}^{(r)}\doteq 1+\sum _{j\in J}{\int }_{{\Gamma }_j}\mathrm{d}{\lambda }_j(k){\varphi }_{n,m}^{(r)}(k){\sigma }_{n,m}^{(r)}({\omega }^{j}k){\omega }^{-j}{k}^{-1},$$

$$p\left({\tau }_{n,m}^{(r+1+\alpha +\beta )}{\tau }_{n+1,m+1}^{(r)}-{\tau }_{n,m+1}^{(r+\alpha )}{\tau }_{n+1,m}^{(r+1+\beta )}\right)=q\left({\tau }_{n,m}^{(r+1+\alpha +\beta )}{\tau }_{n+1,m+1}^{(r)}-{\tau }_{n+1,m}^{(r+\beta )}{\tau }_{n,m+1}^{(r+1+\alpha )}\right),$$

$$p+u_{n,m}^{(r+1)}-u_{n+1,m}^{(r)}=p\frac{{\tau }_{n,m}^{(r)}{\tau }_{n+1,m}^{(r+1)}}{{\tau }_{n,m}^{(r+1)}{\tau }_{n+1,m}^{(r)}},q+u_{n,m}^{(r+1)}-u_{n,m+1}^{(r)}=q\frac{{\tau }_{n,m}^{(r)}{\tau }_{n,m+1}^{(r+1)}}{{\tau }_{n,m}^{(r+1)}{\tau }_{n,m+1}^{(r)}}and{v}_{n,m}^{(r)}=\frac{{\tau }_{n,m}^{(r+1)}}{{\tau }_{n,m}^{(r)}}.$$

$${\tau }_{n,m}^{(r)}=\sum _{i^{\prime }=0}^{\infty }\frac{1}{i^{\prime }!}{\left(\sum _{i=0}^{\infty }\frac{(-1{)}^i}{i+1}\mathrm{tr}{\left(\mathbf{\Omega }{\mathbf{C}}_{n,m}^{(r)}\right)}^{i+1}\right)}^{i^{\prime }},$$

$$\Lambda (i,j)={\delta }_{i+1,j},{}^t\Lambda (i,j)={\delta }_{i,j+1}andO(i,j)={\delta }_{i,0}{\delta }_{0,j};$$

$${\mathbf{C}}_{n,m}^{(r)}=\sum _{j\in J}{\int }_{{\Gamma }_j}\mathrm{d}{\lambda }_j(k){\rho }_{n,m}^{(r)}(k)\mathbf{c}(k){}^t\mathbf{c}({\omega }^{j}k){\sigma }_{n,m}^{(r)}({\omega }^{j}k)$$

$$\mathbf{c}(k)=(\cdots ,k^{-2},k^{-1},1,k,k^2,\cdots {)}^{\mathrm{T}},{}^t\mathbf{c}(k^{\prime })=(\cdots ,k^{\prime -2},k^{\prime -1},1,k^{\prime },k^{\prime 2},\cdots ).$$

Full text

1. Introduction

The study of integrable difference equations is nowadays often considered as a core topic in the theory of integrable systems. Technically, an integrable discrete equation may arise as the superposition formula for Bäcklund transforms of its corresponding continuous equation, see e.g. [23]. However, discrete systems in their own rights have remarkable mathematical significance due to the very rich algebraic structure behind them. The key characteristic of an integrable partial differential equation (PDE) is a Lie algebra, which completely determines the algebraic structure of the equation, see e.g. [7, 14]. Compared with the continuous theory, the integrability of a partial difference equation (P$\Delta$E) is, instead, governed by its associated Lie group. It is such a group structure that makes discrete models algebraically much more general than continuous ones. Reflecting on the nonlinear difference equation itself, this leads to the phenomenon that a suitable continuum limit of a P$\Delta$E obtained from the Bäcklund transform (BT) approach normally yields the whole hierarchy of its corresponding PDEs [24]. From this viewpoint, one could argue that discrete equations are the master models among integrable systems, and continuous equations in a sense are the ‘approximation’ of them. However, most of the mathematical methods for integrable systems are invented based on the continuous theory, and some times are no longer valid to the study of discrete models. Therefore, to completely understand the structure of discrete systems, it is essential to develop classical methods and invent novel concepts for the discrete theory, see monograph [13].

A milestone in the study of discrete integrable systems is the Adler–Bobenko–Suris (ABS) classification [1] of scalar integrable quadrilateral equations (i.e. P$\Delta$Es defined on four points). In this classification, the most well-known model is the discrete Korteweg–de Vries (KdV) equation (see Figure 1(a))

[TABLE]

Any of the above $\tau$-equations can be considered as the discrete bilinear KdV equation, since not only the continuum limit of each equation gives rise to the continuous bilinear KdV equation, but also the tau functions are exactly the same from the viewpoint of soliton-type solution, describing structure of the discrete KdV equation, see [13]. However, we still needs to mark that from the aspect of initial value problem of lattice equations, the 5-point equation and the 6-point equations are not equivalent.

Apart from quadrilateral equations, there also exist a large class of integrable difference equations defined on a $3\times 3$ stencil (i.e. 9-point equations which are second order in both $n$ and $m$), known as the discrete Boussinesq (BSQ) family. A typical model is the discrete BSQ equation (see Figure 2)

[TABLE]

This may look surprising since the continuous BSQ equation can be written as a bilinear equation. However, if we take a suitable continuum limit of the 9-point equation (1.2), a common factor is generated in the procedure, and thus, the continuum limit of the trilinear equation does give rise to the continuous bilinear BSQ equation. This gives us a strong hint that in the discrete theory, the bilinear structure may no longer be preserved. The discrete BSQ-type equations can alternatively be written as three-component (effectively two-component) coupled systems of first-order P$\Delta$Es. A search for integrable discrete BSQ-type systems was made by Hietarinta [12], resulting in a remarkable classification of integrable multi-component BSQ-type systems. It was later proven that all Hietarinta’s discrete BSQ-type systems arise from the so-called extended discrete BSQ systems [26].

The quadrilateral KdV-type equations and 9-point BSQ systems are integrable in the sense of the so-called multi-dimensional consistency (MDC) property (i.e. a P$\Delta$E can be consistently embedded into a higher-dimensional lattice space, see e.g. Doliwa and Santini [6], Nijhoff and Walker [21], and also Bobenko and Suris [4]). The MDC of the integrable quadrilateral equations is often presented as the consistency-around-the-cube (CAC) phenomenon, namely embedding a quadrilateral equation into a cube and showing the compatibility of equations defined on different faces, see e.g. [13]. However, the verification of the CAC property is sometimes subtle. The famous discrete sine–Gordon (sG) equation is a typical example – one must consider simultaneously the discrete sG and the discrete modified KdV on a cube, rather than verify its consistency directly; to be more precise, we have the discrete sG on the top and bottom faces and the discrete modified KdV on the four side faces on a single cube. In other words, the discrete mKdV and sG flows commute. Therefore, the key point behind the MDC property is actually the existence of infinitely many commuting discrete flows of a difference equation, which in turn results in the integrability (e.g. existence of a Lax pair, see [17, 4]) of a difference equation. More generally, an integrable discrete equation is compatible with not only its discrete symmetries, but also its semi-discrete (semi-discrete hierarchy) and continuous symmetries (i.e. the continuous hierarchy); in other words, discrete and continuous flow variables are on the same footing, and the consistency is guaranteed by a unified solution structure behind all the equations, see [11].

The discrete KdV-type and BSQ-type equations are lower-order lattice models. A natural problem would be to understand the integrability of higher-order (or higher-rank) P$\Delta$Es. However, studying higher-order difference equations directly is rather difficult since even the 9-point P$\Delta$Es already take complex forms. A sensible way to solve such a problem is, instead, to consider integrable discrete multi-component systems, which, in a sense, makes it possible for us to deal we equations in relatively lower orders, though more components must be introduced as the price we have to pay. Typical ways of constructing such systems include periodic reductions of discrete higher-dimensional models (see e.g. [15]) and the Lax approach based on “big” matrices (for instance, see [22]), etc.

Very recently, Fordy and Xenitidis [9] generalised the early result regarding the discrete-time Bogoyavlensky lattice equations in [22], and proposed a systematical scheme for integrable difference equations associated with $\mathbb{Z}_{\mathcal{N}}$ graded discrete Lax pairs. The key point in their construction is an $\mathcal{N}\times\mathcal{N}$ periodic matrix

[TABLE]

which has very nice properties such as $\mathbf{\Sigma}^{\mathcal{N}}=\mathbf{I}$ and $\mathbf{\Sigma}^{\mathrm{T}}=\mathbf{\Sigma}^{-1}$, where $\mathbf{I}$ is the $\mathcal{N}\times\mathcal{N}$ unit matrix. With the help of $\mathbf{\Sigma}$, $\mathbb{Z}_{\mathcal{N}}$ graded discrete Lax pairs are generated by the $\mathcal{N}\times\mathcal{N}$ level $\alpha$ matrices $\mathbf{A}=\operatorname{diag}(a^{(0)},a^{(1)},\cdots,a^{(\mathcal{N}-1)})\mathbf{\Sigma}^{\alpha}$ for $\alpha=0,1,\cdots,\mathcal{N}-1$, i.e. the $\mathbb{Z}_{\mathcal{N}}$ graded ring. The benefit of such a construction is that all the integrable difference equations arising from such a framework can be written as quadrilateral coupled systems and one thus successfully avoids dealing with higher-order terms. Moreover, the integrability of these integrable difference equations also results in novel Yang–Baxter (YB) maps [10].

However, the exact solutions to the Fordy–Xenitidis (FX) models still needs to be clarified. Our motivation is to investigate solution structure of these equations. In this paper, we consider discrete equations in the equivalence class $[(\alpha,\alpha+1);(\beta,\beta+1)]$ within the FX classification. Integrable discrete equations in this case can mainly be expressed by either of the $\mathcal{N}$-component discrete systems111 It is obvious that in equation (1.4) the lattice parameters $p$ and $q$ can be absorbed by introducing a simple point transform (i.e. a translation), but here we still reserve the lattice parameters in the equation. This is because this treatment allows us to establish a more direct one-to-one correspondence between the system and its continuous analogue, and as a by-product, the procedure of taking continuum limit is also more natural.

[TABLE]

and

[TABLE]

in which $r=0,1,\cdots,\mathcal{N}-1$ for $\mathcal{N}=2,3,\cdots$. The variables $u$ and $v$ are referred to as the additive and quotient potentials, respectively, describing the structure of unmodified and modified equations, respectively. It was also pointed out in reference [9] that the additive potential $u$ and the quotient potential $v$ satisfy identities

[TABLE]

respectively.

These identities theoretically reduce the number of components in the discrete systems (1.4) and (1.5). For instance, when $(\alpha,\beta)=(0,0)$, $u_{n,m}^{(\mathcal{N}-1)}$ and $v_{n,m}^{(\mathcal{N}-1)}$ can be replaced by the other components through these identities, which result in $(\mathcal{N}-1)$-component systems

[TABLE]

for the additive potential and

[TABLE]

for the quotient potential, i.e. the discrete GD hierarchy (see [20]) and the discrete modified GD hierarchy (see e.g. [3, 5]); similarly, in the $(\alpha,\beta)=(\mathcal{N}-1,\mathcal{N}-1)$ case, we have

[TABLE]

namely the discrete Schwarzian GD hierarchy which was first proposed in [2]. When $\mathcal{N}=2$, the above three systems can be further rewritten as scalar equations (1.1b), (1.1a) and (1.1c), respectively; while for $\mathcal{N}=3$, they lead to the 9-point equations (1.2), (1.2) and (1.2), respectively.

Formulae (1.4) and (1.5) together with (1.6) not only cover the well-known discrete models mentioned above, but also give rise to a huge number of new integrable multi-component discrete systems in case of $\mathcal{N}\geq 3$, see [9] for more details. For example, when $\mathcal{N}=3$ and $(\alpha,\beta)=(1,1)$, the multi-component system of the quotient potential $v=v^{(0)}$ can reduce to a scalar integrable discrete equation

[TABLE]

where $X$ and $Y$ are given by

[TABLE]

respectively, which is a novel 9-point lattice equation compared with (1.2), (1.2) and (1.2).

Both systems (1.4) and (1.5) are integrable in the sense of having their respective Lax pair. Concretely, the compatibility condition of the $\mathcal{N}$-component linear system (i.e. $\mathbb{Z}_{\mathcal{N}}$ graded Lax pair)

[TABLE]

in which $\phi_{n,m}^{(r+\mathcal{N})}=k^{\mathcal{N}}\phi_{n,m}^{(r)}$ and $r=0,1,\cdots,\mathcal{N}-1$ for $\mathcal{N}=2,3,\cdots$, gives rise to equation (1.4). Furthermore, the potentials $u$ and $v$ are connected via a difference transform [9], i.e. a BT, which is composed of

[TABLE]

Substituting $u$ with $v$ with the help of (1.9) in (1.8), we obtain the Lax pair for equation (1.5).

The method we adopt is the so-called direct linearisation (DL), which was originally proposed by Fokas and Ablowitz in [8] to solve the initial value problem of the KdV equation. It was later developed by Nijhoff, Quispel, Capel et al. to construct integrable discrete equations and their exact solutions, see e.g. [19, 20]. The advantage of the DL is that it not only provides solutions to nonlinear equations, but also reveals the underlying structures of nonlinear integrable equations. The main idea in the DL approach is to associate a nonlinear equation with a linear integral equation. In practice, by considering the nonlinearisation of a certain linear integral equation, one can algebraically construct the corresponding nonlinear equation as well as its solutions and other integrability characteristics simultaneously.

The main results in this paper are listed below.

Theorem 1**.**

Consider a linear integral equation

[TABLE]

where $\rho$ and $\sigma$ are the plane wave factors are expressions of the discrete variables $n$, $m$ and $r$ (which are associated with their respective lattice parameters $p$, $q$ and [math]) as well as the spectral parameter $k$ (or $l$), given by

[TABLE]

$\phi$* is the wave function depending on the discrete variables $n$, $m$, $r$ and also the spectral parameter; $\mathrm{d}\lambda_{j}$ are the measures for the integrations on their respective contours $\Gamma_{j}$ (which are arbitrary and can be specified for particular classes of solutions); $\omega\doteq\exp\{2\pi\mathfrak{i}/\mathcal{N}\}$ with $\mathfrak{i}$ being the imaginary unit is a primitive $\mathcal{N}$th root of unity; and $J$ is a set defined as J\doteq\{j|\hbox{0<j<\mathcal{N}$are integers coprime to$\mathcal{N}}\}.*

For any solution $\phi$ to the linear integral equation (1.10), the additive and quotient potentials given by

[TABLE]

respectively, solve equations (1.4) and (1.5), respectively, and $\phi$ solves the corresponding linear system (1.8).

Theorem 2**.**

Equations (1.4) and (1.5) are both linked with the bilinear system

[TABLE]

in which $\tau_{n,m}^{(r+\mathcal{N})}=\tau_{n,m}^{(r)}$, and $r=0,1,\cdots,\mathcal{N}-1$ for $\mathcal{N}=2,3,\cdots$ through

[TABLE]

The bilinear system (1.13) is solved by

[TABLE]

where $\mathbf{\Omega}$ is an infinite matrix defined as $\mathbf{\Omega}=\sum_{i=0}^{\infty}{}^{t\!}\mathbf{\Lambda}^{i}\mathbf{O}\mathbf{\Lambda}^{-i-1}$ with $\mathbf{\Lambda}$, ${}^{t\!}\mathbf{\Lambda}$ and $\mathbf{O}$ being infinite matrices having their respective $(i,j)$-entries

[TABLE]

the infinite matrix $\mathbf{C}$ is defined as

[TABLE]

in which $\mathbf{c}$ and ${}^{t\!}\mathbf{c}$ are an infinite-dimensional column vector and an infinite-dimensional row vector, respectively, defined as follows:

[TABLE]

Theorem 1 establishes the link between equations (1.4) and (1.5), and the linear integral equation (1.10). This theoretically provides us with a way to solve the two discrete nonlinear systems. Theorem 2 gives us not only a novel nonlinear form (i.e. the bilinear form) of the FX equations, but also its solution.

2. Proof of theorem 1

We consider a more general linear integral equation

[TABLE]

in which the wave function $\Phi$ is an infinite column vector with $i$th component $\phi_{n,m}^{(h)}(i;k)$ as a function depending on the discrete variables $n$, $m$ and $h$ as well as the spectral variable $k$ (or $l$); and the plane wave factors are given by

[TABLE]

namely the $(\alpha,\beta)=(0,0)$ case of (1.11). Compared with (1.10), (2.1) is a set of linear integral equations. The reason why we start with such a linear integral equation is because it allows us to introduce the so-called infinite matrix formalism, which will significantly simplify the algebraic construction of discrete nonlinear equations. The idea of the whole proof is to link equation (2.1) with (1.4) and (1.5) for the $(\alpha,\beta)=(0,0)$, and then by introducing the parameter deformation with respect to $\alpha$ and $\beta$, we complete the proof of theorem 1.

We introduce a new quantity

[TABLE]

The quantity $\mathbf{U}$ is an infinite matrix having its $(i,j)$-entry $U_{n,m}^{(h)}(i,j)$ only depending on the discrete variables $n$, $m$ and $h$, and independent of the spectral parameter. It plays the role of nonlinearisation of the wave function $\Phi$. On the other hand, following the definition of $\mathbf{\Omega}$, we can verify that

[TABLE]

directly. Replacing $\frac{1}{k-\omega^{j}l}$ by ${}^{t\!}\mathbf{c}(\omega^{j}l)\mathbf{\Omega}\mathbf{c}(k)$ in the integral equation (2.1), we have

[TABLE]

in other words, equation (2.1) can be reformulated as a compact form

[TABLE]

Taking the definition of $\mathbf{U}$ (i.e. (2.3)) into consideration and acting the operation $\sum_{j\in J}\int_{\Gamma_{j}}\mathrm{d}\lambda_{j}(k)\cdot{}^{t\!}\mathbf{c}(\omega^{j}k)\sigma_{n,m}^{(h)}(\omega^{j}k)$ on the linear integral equation (2.4), we immediately obtain its nonlinear version, i.e. a relation for the quantity $\mathbf{U}$, as follows:

[TABLE]

where $\mathbf{C}$ is an infinite matrix defined as

[TABLE]

which is actually the infinite matrix representation of the effective plane wave factor

[TABLE]

We now algebraically construct equations (1.4) and (1.5) based on the infinite matrix $\mathbf{U}$. Since $\omega^{j}$ has the property $(\omega^{j})^{\mathcal{N}}=(\omega^{\mathcal{N}})^{j}=1$, we can easily prove from (2.7) that the effective plane wave factor satisfies the periodic condition $\rho_{n,m}^{(h+\mathcal{N})}(k)\sigma_{n,m}^{(h+\mathcal{N})}(\omega^{j}k)=\rho_{n,m}^{(h)}(k)\sigma_{n,m}^{(h)}(\omega^{j}k)$. This immediately results in the fact that the infinite matrix $\mathbf{C}$ obeys the periodic condition

[TABLE]

as the dynamics of $\mathbf{U}$ relies on $\mathbf{C}$ explicitly according to (2.5), which implies that the discrete variable $h$ from now on can only be selected from a finite set, namely it will play the role of labelling components in the resulting nonlinear equations. Notice that

[TABLE]

as well as $\mathbf{\Omega}\mathbf{\Lambda}-{}^{t\!}\mathbf{\Lambda}\mathbf{\Omega}=\mathbf{O}$. By shifting equation (2.5) with respect to $n$ by one unit, one obtains

[TABLE]

which can equivalently be written as

[TABLE]

Multiplying this equation by $\left(1+\mathbf{\Omega}\mathbf{C}_{n,m}^{(h)}\right)^{-1}$ and noticing equation (2.5), we finally reach

[TABLE]

Equations listed in (2.9) are the fundamental relations to algebraically construct the nonlinear equations. The main idea of the construction is to focus on certain entries of $\mathbf{U}$ and prove that these entries form their respective closed-form nonlinear equations, which will later lead to equations (1.4) and (1.5). For this purpose, we introduce new variables

[TABLE]

It is obvious to see that $u$ and $v$ satisfy

[TABLE]

because they are effectively the entries of $\mathbf{U}$ and the periodicity is preserved. Taking the $(0,0)$-entry of (2.9), the dynamical equations in (2.9) give rise to the following equations involving only $u$, $U(1,0)$ and $U(0,1)$:

[TABLE]

We aim to eliminate $U(1,0)$ and $U(0,1)$ and derive a closed-form equation for only $u$. In practice, by subtraction the entry $U(1,0)$ can be eliminated and we obtain

[TABLE]

The summation $\eqref{dGD:a}_{n,m}^{(h+1)}+\eqref{dGD:b}_{n+1,m}^{(h)}+\eqref{dGD:c}_{n,m+1}^{(h)}$ can help to eliminate the entry $U(0,1)$, and thus, we end up with a closed-form discrete equation expressed by only $u$, which is in the form of

[TABLE]

This equation together with the periodic condition for $u$ given in (2.11), namely $u_{n,m}^{(h+\mathcal{N})}=u_{n,m}^{(h)}$, forms the multi-component system (1.4) for $(\alpha,\beta)=(0,0)$, when $h$ is fixed from $0,1,\cdots,\mathcal{N}-1$. Next, we construct the closed-form equation for the potential $v$. The $(0,-1)$-entry of (2.9) gives us the following:

[TABLE]

Eliminating $U(1,-1)$ in these relations, the following difference equations involving $u$ and $v$ can be derived:

[TABLE]

which actually play the role of discrete Miura-type transforms, though the transforms here are nonlocal in both directions. Adding the three equations in (2.14) up, the additive potential $u$ is eliminated, and it turns out that we derive a closed-form equation for only $v$, which is in the form of

[TABLE]

Notice that $v$ satisfies the periodicity (2.11). Equation (2.15) together with $v_{n,m}^{(h+\mathcal{N})}=v_{n,m}^{(h)}$ for $h=0,1,\cdots,\mathcal{N}-1$ forms the $(\alpha,\beta)=(0,0)$ case of (1.5).

Next, we prove that the component

[TABLE]

is the wave function for the Lax pairs of equations (2.13) and (2.15), respectively. This can be realised by constructing the Lax pairs by starting from (2.4). Remind ourselves that $\mathbf{U}$ obeys the periodic condition $\mathbf{U}_{n,m}^{(h+\mathcal{N})}=\mathbf{U}_{n,m}^{(h)}$ and the plane wave factor $\rho$ satisfies $\rho_{n,m}^{(h+\mathcal{N})}(k)=k^{\mathcal{N}}\rho_{n,m}^{(h)}(k)$. The relation (2.4) yields the quasi-periodic condition

[TABLE]

because $\phi$ is a component of $\Phi$. Shifting the wave function $\Phi$ with respect to $n$ by one unit, we have

[TABLE]

which implies that the wave function $\Phi$ satisfies a dynamical relation

[TABLE]

Subtracting (2.17c) from (2.17a) and (2.17b), respectively, we are able to eliminate $\mathbf{\Lambda}$ and obtain the following dynamical relations:

[TABLE]

These two equations will lead to the Lax pair for equation (2.13). In fact, by taking the [math]-component of the two equations, we immediately get

[TABLE]

which together with the quasi-periodic condition $\phi_{n,m}^{(h+\mathcal{N})}=k^{\mathcal{N}}\phi_{n,m}^{(h)}$ forms the $\mathbb{Z}_{\mathcal{N}}$ graded Lax pair for (2.13), for $h=0,1,\cdots,\mathcal{N}-1$. Replacing $u$ by $v$ with the help of the discrete Miura transforms given in (2.14), the Lax pair for equation (2.15), namely

[TABLE]

in which $h=0,1,\cdots,\mathcal{N}-1$, is obtained.

We have proven the $(\alpha,\beta)=0$ case of theorem 1, and now consider the case for generic $\alpha$ and $\beta$. Comparing the plane wave factors (2.2) and (1.11), we can observe that the latter is a parameter deformation of the former, by introducing a change of variables $(n,m,h)\mapsto(n,m,r)$ as follows:

[TABLE]

where $\alpha,\beta\in\{0,1,\cdots,\mathcal{N}-1\}$. The parameter deformation results in the dynamical relations $\rho_{n+1,m}^{(r)}=\rho_{n+1,m}^{(h+\alpha)}$, $\rho_{n,m+1}^{(r)}=\rho_{n,m+1}^{(h+\beta)}$, and $\rho_{n,m}^{(r+1)}=\rho_{n,m}^{(h+1)}$ for $\rho$ as well as the same relations for $\sigma$. Since the dynamics of $\mathbf{C}$ explicitly relies on the product of $\rho$ and $\sigma$, one can verify that the infinite matrix $\mathbf{C}$ satisfies

[TABLE]

as they follow from (2.6) and (1.16) directly. Furthermore, the infinite matrix $\mathbf{U}$ obeys relations $\mathbf{U}_{n+1,m}^{(r)}=\mathbf{U}_{n+1,m}^{(h+\alpha)}$, $\mathbf{U}_{n,m+1}^{(r)}=\mathbf{U}_{n,m+1}^{(h+\beta)}$ and $\mathbf{U}_{n,m}^{(r+1)}=\mathbf{U}_{n,m}^{(h+1)}$, as the dynamics of $\mathbf{U}$ only depends on $\mathbf{C}$ through the equation $\mathbf{U}=(1-\mathbf{U}\mathbf{\Omega})\mathbf{C}$, cf. (2.5), which consequently together with the $\rho$ relations result in the same deformation on the wave function $\Phi$, namely $\Phi_{n+1,m}^{(r)}=\Phi_{n+1,m}^{(h+\alpha)}$, $\Phi_{n,m+1}^{(r)}=\Phi_{n,m+1}^{(h+\beta)}$, and $\Phi_{n,m}^{(r+1)}=\Phi_{n,m}^{(h+1)}$. Noticing that the additive potential $u$ and the quotient potential $v$, and the wave function $\phi$ are entries and the [math]th-component of $\mathbf{U}$ and $\Phi$, respectively, see (2.10) and (2.16), we conclude that these variables obey exactly the same deformation relations. And thus, by applying the deformation to nonlinear equations (2.13) and (2.15), the difference transform (2.14), and the Lax pair (2.19), one obtains the FX discrete models (1.4) and (1.5), the relevant BT (1.9), as well as the $\mathbb{Z}_{\mathcal{N}}$ graded Lax pair (1.8). To put it in another way, based on the linear integral equation

[TABLE]

we have proven that $u=U(0,0)$ and $v=1+U(0,-1)$ solve the nonlinear equations (1.4) and (1.5), and $\phi(k)=\Phi(0,k)$ solves the associated Lax pair (1.8), which is indeed the result of theorem 1.

3. Proof of theorem 2

We still start with the $(\alpha,\beta)=(0,0)$ case. The idea is to establish the link between dynamics of the tau functions and entries of the infinite matrix $\mathbf{U}$. By making use of the infinite matrix relations found in section 2, a closed-form multicomponent system expressed by only the tau functions, namely the bilinear form of the FX nonlinear discrete systems is therefore obtained.

The tau function in this case is defined as $\tau_{n,m}^{(h)}=\det\left(1+\mathbf{\Omega}\mathbf{C}_{n,m}^{(h)}\right)$, in which the determinant of infinite matrix should be thought of as the expansion of

[TABLE]

as the determinant of an infinite matrix is defined thorough $\ln\det(\mathbf{A})=\operatorname{tr}\ln(\mathbf{A})$. A comment here is that the infinite summation does not lead to divergence, because the infinite matrix $\mathbf{\Omega}$ depends on $\mathbf{O}$ and the trace operation guarantees the convergence. From the definition of the tau function, one can immediately deduce that the tau function satisfies the periodic condition

[TABLE]

This is because the dynamics of $\tau$ is completely governed by $\mathbf{C}$ according the definition (3.1), and therefore the tau function inherits the periodicity of $\mathbf{C}$, cf. (2.8).

Next, we consider the dynamics of the tau function. For future convenience, we introduce a new quantity

[TABLE]

where $(\cdot)(0,0)$ denotes taking the $(0,0)$-entry of the infinite matrix in the bracket, and $a$ is a parameter. Now we establish the link between the dynamics of the tau function and the entries of the infinite matrix $\mathbf{U}$. Direct calculation shows that

[TABLE]

where the Weinstein–Aronszajn formula is used for the second equality. This gives rise to

[TABLE]

These equalities form the fundamental relations to construct the bilinear form of the FX discrete equations. Multiplying equation (2.9a) by $(a+{}^{t\!}\mathbf{\Lambda})^{-1}$, the following equation is obtained:

[TABLE]

Taking the $(0,0)$-entry of this matrix equation, we can derive

[TABLE]

Similarly, we have also

[TABLE]

for the $m$-direction. Subtracting the above two equations, the entry $\left(\mathbf{\Lambda}\mathbf{U}(a+{}^{t\!}\mathbf{\Lambda})^{-1}\right)(0,0)$ is eliminated, and thus, we end up with a difference relation between $u$ and $V(a)$, which is given by

[TABLE]

By setting $a=q$ and taking the first two equations in (3.3) into consideration, we can replace the $V$-variable by the tau function, and consequently obtain

[TABLE]

Eliminating $u$ in (3.4), a bilinear equation arises, which takes the form of

[TABLE]

Equation (3.5) subject to the periodic relation (3.2) for $h=0,1,\cdots,\mathcal{N}-1$ forms an $\mathcal{N}$-component coupled discrete bilinear system, which is the bilinear version of the FX models when $(\alpha,\beta)=(0,0)$.

We can now apply the parameter deformation to the bilinear system (3.5). Notice that in (1.15) the dynamics of the tau function $\tau$ is completely determined by $\mathbf{C}$, we thus have a similar change of variables, namely

[TABLE]

The application of equation (3.6) to (3.5) gives us the bilinear equation (1.13). Simultaneously, the tau function (3.1) turns out to be (1.15). Moreover, we note that the acting the same change of variables on (3.4) and (3.3) also provides us with (1.14), which together with the periodicity condition of the tau function, namely $\tau_{n,m}^{(r+\mathcal{N})}=\tau_{n,m}^{(r)}$, yield the constraints (1.6). Theorem 2 is proven.

4. Concluding remarks

In this paper, theorems 1 and 2 provide solutions involving an integral with respect to the spectral variable to the FX models (1.4) and (1.5), and the bilinear equation (1.13), respectively. Our main motivation is to explain the solution structure of the integrable difference equations proposed by Fordy and Xenitidis. In fact, in the solution structure, the spectral variable sits in a continuous space; by specifying the measure associated with the linear integral, special classes of solutions to these discrete systems will naturally arise. Particularly, reduction of the measure to a discrete space, i.e. introducing a finite number of singularities in the measure, will lead to the degeneracy of the linear integral, providing determinantal (namely soliton-type) solutions.

We take equation (1.7) as an example. In this case $\mathcal{N}=3$, $\alpha=\beta=1$, and the linear integral equation (1.10) becomes

[TABLE]

where the plane wave factors satisfy

[TABLE]

and $\omega=\exp\{2\pi\mathfrak{i}/3\}$. Accordingly, by following (1.12), we have

[TABLE]

Now we take very special measures

[TABLE]

and require that the contours $\Gamma_{1}$ and $\Gamma_{2}$ contain $k_{1}$ and $k_{2}$ as their only singularities, respectively. As a result, by the residue theorem equation (4.1) yields

[TABLE]

which in turn implies

[TABLE]

if $k$ is set to be $k_{1}$ and $k_{2}$, respectively. Simultaneously, the special measure also give rise to a reduction of the potential $v$ given by

[TABLE]

And thus, we obtain

[TABLE]

Although the above formula contain two parameters $k_{1}$ and $k_{2}$, we note that for any $r=0,1,2$ the obtained $v$ plays the role of the one-soliton solution to the FX equation. This is because the two parameters really follow from the two distinct measures and they do not govern soliton interaction; in other words, the result obtained from the DL method in a sense generalises the traditional one-soliton solution.

Another remark is that for fixed $\mathcal{N}$, all FX discrete equations share the same form of linear integral equation (1.10) but with different plane wave factors (1.11). In other words, from the viewpoint of solution structure, the FX systems are distinguished by the linear dispersion, while the nonlinear structures (i.e. the behaviours of the spectral variable in difference equations) are the same. As the $(\alpha,\beta)=(0,0)$ case describes the structure of the standard discrete GD hierarchy, the novel equations in the FX in a sense can also be thought of as the equations governed by deformed discrete flows in the discrete GD hierarchy.

We only consider the coprime case of the FX classification in this case. It is not yet clear whether discrete equations in the non-coprime case can also be solved from the DL approach, but it is quite certain that solution structure of those equations is entirely different from that of the discrete GD hierarchy at least. This remains a problem for our future work.

Acknowledgments

This project was sponsored by the National Natural Science Foundation of China (grant no. 11901198) and by Shanghai Pujiang Program (grant no. 19PJ1403200). The author was also partially supported by the Science and Technology Commission of Shanghai Municipality (grant no. 18dz2271000), and he thanks Frank Nijhoff for comments on a draft version of this manuscript.