# Integrability conditions for Boussinesq type systems

R. Hernández Heredero  ETSIS de Telecomunicación, Universidad Politécnica de
Madrid, Madrid, Spain  and  V. Sokolov  Higher School of Modern Mathematics
MIPT, Moscow, Russia

(Date: )

###### Abstract.

The symmetry approach to the classification of evolution integrable partial
differential equations (see, for example [2]) produces an infinite series of
functions, defined in terms of the right hand side, that are conserved
densities of any equation having infinitely many infinitesimal symmetries. For
instance, the function $\frac{\partial f}{\partial u_{x}}$ has to be a
conserved density of any integrable equation of the KdV type
$u_{t}=u_{xxx}+f(u,u_{x})$. This fact imposes very strong conditions on the
form of the function $f$. In this paper we construct similar canonical
densities for equations of the Boussinesq type. In order to do that, we write
the equations as evolution systems and generalise the formal diagonalisation
procedure proposed in [1] to these systems.

##  1\. Introduction

Integrable systems constitute an important class of equations because they
have multiple applications in mathematical physics and pure mathematics. They
are present in fluid mechanics, nonlinear optics, field theory, numerical
analysis, algebraic geometry, etc., and explain phenomena such as the soliton
and instantons, besides appearing in spectral theory and the theory of
differential and pseudo-differential operators. Its study boomed since the
discovery of the inverse scattering transform in the 1960’ies, and has not
stopped since then.

The main feature of integrable systems is that they possess an infinite number
of non-trivial symmetries and a bihamiltonian structure that, in nonlinear
systems, relates the symmetries to an infinite number of non-trivial
conservation laws.

The _symmetry approach_ (cf. [2, 3]) provides a powerful, algorithmic method
to study integrability, defining ways to derive important objects such as
recursion, symplectic or Hamiltonian operators admitted by integrable systems.
The symmetry approach has been extensively applied to scalar evolution
equations, such as those of KdV type, systems of evolution equations like NLS
type equations, or Boussinesq type systems of equations. As the complexity of
the system under study grows, for example in the number of independent or
dependent variables involved, or for non-evolutionary equations, the approach
has to be further developed. In this paper we have systematised and extended
the symmetry approach to the study of equations of the Boussinesq type,
defining the concept of regularly diagonalisable systems and providing a
general formula for canonical conserved densities of certain subclass of
Boussinesq type of systems. We write explicitly the first canonical
conservation laws of a subclass of fourth order Boussinesq type family of
equations, that we finally classify.

Non-evolution (1+1)-dimensional PDEs of the form

(1) $u_{tt}=f(u,u_{1},\ldots,u_{n};u_{t},u_{1t},\ldots,u_{mt}),\qquad
u_{i}=\frac{\partial^{i}u}{\partial x^{i}},\quad
u_{it}=\frac{\partial^{i+1}u}{\partial x^{i}\partial t}$

are called equations of Boussinesq type. Such equations can be rewritten as
evolution systems of two variables as

(2) $\displaystyle u_{t}$ $\displaystyle=v,$ $\displaystyle v_{t}$
$\displaystyle=f(u,u_{1},\ldots,u_{n},v,v_{1},\ldots,v_{m}).$

In this paper we consider the class of integrable111In our paper integrability
means the existence of local infinitesimal symmetries and/or conservation laws
(see, for example, [4, 5, 6]). systems (2) of the form

(3) $\begin{aligned} u_{t}&=v,\\\
v_{t}&=f(u,u_{1},\ldots,u_{2k},v,v_{1},\ldots,v_{k-1}),\end{aligned}\qquad
k>1.$

The Boussinesq equation $u_{tt}=u_{4}+2u_{1}u_{2}$, written as the system

$\displaystyle u_{t}$ $\displaystyle=v,$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+2u_{1}u_{2}$

belongs to this class.

######  Remark 1.

A similar class of integrable systems

(4) $\displaystyle u_{t}$ $\displaystyle=v_{1},$ $\displaystyle v_{t}$
$\displaystyle=u_{2k-1}+g$

was studied in [7]. If $f_{u}=0$, then setting $u_{1}\to\tilde{u}$ allows to
write (3) as

$\displaystyle\tilde{u}_{t}$ $\displaystyle=v_{1},$ $\displaystyle v_{t}$
$\displaystyle=f(\tilde{u},\ldots,\tilde{u}_{2k-1},v,v_{1},\ldots,v_{k-1}).$

If $f=Dg$ (where $Dg$ will denote the total $x$-derivative of a differential
function $g$ throughtout the paper), setting $v\to\tilde{v}_{1}$ leads to

$\displaystyle u_{t}$ $\displaystyle=\tilde{v}_{1},$
$\displaystyle\tilde{v}_{t}$
$\displaystyle=g(u,u_{1},\ldots,u_{2k-1},\tilde{v}_{1},\ldots,\tilde{v}_{k})$

so the corresponding function $g$ on (4) must not depend on $v$. Thus both
families of systems (3) and (4) have an intersection but are not equivalent.

In the paper [8] integrable systems of the form

$\displaystyle u_{t}$ $\displaystyle=v_{1},$ $\displaystyle v_{t}$
$\displaystyle=\alpha u_{2k-1}+\beta
v_{k}+g(u,\ldots,u_{2k-2},v,\ldots,v_{k-1})$

with constant $\alpha$ and $\beta$ were considered. A polynomiality
restriction on $g$ was crucial due to the symbolic approach used in that
paper. A diagonalisation of the linear part (separant) of the system, in
symbolic representation, was used to perform a complete classification of
homogeneous polynomial systems of that type with $k=2,3,5$. Some non-
polynomial examples of equations (1) with $k=3$ were found in [9].

The symmetry approach [2, 10, 11] to the classification of integrable scalar
evolution equations

(5) $u_{t}=\phi(u,u_{1},\dots,u_{n}),\qquad
u_{i}=\frac{\partial^{i}u}{\partial x^{i}}$

is based on the existence of a formal pseudo-differential recursion operator
of the form

(6) $R=r_{1}D+r_{0}+r_{-1}D^{-1}+\cdots,\qquad r_{i}=r_{i}(u,u_{1},\dots)$

satisfying, by definition, the following operator relation

(7) $R_{t}=[\phi_{*},\,R].$

Here

$\phi_{*}=\sum_{i=0}^{n}\frac{\partial\phi}{\partial u_{i}}D^{i}$

is the Fréchet derivative of $\phi$, $D$ is the total $x$-derivative, and the
$t$-derivative in (7) is found by virtue of (5). A formal recursion operator
(6) produces [11] an infinite sequence of so-called local canonical conserved
densities

(8) $\rho_{i}={\rm res}\,R^{i},\quad i=-1,1,2,\ldots,\qquad\rho_{0}={\rm
res}\,\log(R).$

Requiring locality of the canonical conserved densities leads to efficient
necessary integrability conditions upon the right hand side $\phi$ (cf. [12]).

Equations (5) that have an infinite sequence of local symmetries of the form

(9) $u_{\tau}=\psi(u,u_{1},\dots,u_{i})$

possess a formal recursion operator [10, 11].

If the equation has local conservation laws then, besides formal recursion
operators, it admits [13] a formal symplectic operator $S$ such that

(10) $S_{t}+S\,\phi_{*}+\phi_{*}^{+}\,S=0\,,$

where + denotes formal operator conjugation. The existence of a formal
symplectic operator $S$ implies that the density $\rho_{2i}$ is trivial for
any $i$.

The series $H=S^{-1}$ satisfies the relation

(11) $H_{t}-H\phi_{*}^{+}-\phi_{*}\,H=0\,,$

and is called formal Hamiltonian operator.

######  Remark 2.

Suppose that both formal recursion and symplectic operators $R$ and $S$ exist.
Since $S^{+}$ and $SR^{i}$ are also formal symplectic operators for any
$i\in\mathbb{Z}$, we may assume [11] that $R$ and $S$ are first order series
such that

(12) $R^{+}=-S\,R\,S^{-1},\qquad{S}^{+}=-S.$

Since ${\rm res}\,X=-{\rm res}\,X^{+}$ and ${\rm res}\,(QXQ^{-1}-X)\in{\rm
Im}\,D$ for any series $X,Q$, it follows from (12) that ${\rm
res}\,R^{2i}\in{\rm Im}\,D$ for any $i$.

######  Definition 1.

Equations (5) that have both formal recursion and symplectic operators are
called S-integrable.

In the case of general multi-component evolution systems

(13) ${\bf u}_{t}=\mathbf{\Phi}({\bf u},{\bf u}_{1},\dots,{\bf u}_{n}),$

where ${\bf u}=(u^{1},\dots,u^{m})$, and
$\mathbf{\Phi}=(\phi^{1},\ldots,\phi^{m})$, the Fréchet derivative
$\mathbf{\Phi}_{*}$ is a matrix of differential operators

(14) $\mathbf{\Phi}_{*}[i,j]=\sum_{l=0}^{n}\frac{\partial\phi^{i}}{\partial
u_{l}^{j}}D^{l}=\mathbf{\Sigma}\,D^{n}+\sum_{l=0}^{n-1}\mathbf{\Phi}^{(l)}D^{l},$

where $\mathbf{\Sigma}$ is the so called _separant_ of the evolution system.
Formal recursion $\mathbf{R}$, symplectic $\mathbf{S}$ and Hamiltonian
$\mathbf{H}$ operators are matrix pseudo-differential series satisfying the
equations

(15) $\displaystyle\mathbf{R}_{t}=[\mathbf{\Phi}_{*},\,\mathbf{R}],$ (16)
$\displaystyle\mathbf{S}_{t}+\mathbf{S}\,\mathbf{\Phi}_{*}+\mathbf{\Phi}_{*}^{+}\,\mathbf{S}=0,$
(17)
$\displaystyle\mathbf{H}_{t}-\mathbf{H}\mathbf{\Phi}_{*}^{+}-\mathbf{\Phi}_{*}\mathbf{H}=0$

where + denotes matrix transposition followed by the formal operator
conjugation of entries.

######  Definition 2.

A system (13) is called non-degenerate if its separant matrix
$\mathbf{\Sigma}$ is invertible and has no multiple eigenvalues at a generic
point.

The symmetry approach was generalised to the case of non-degenerate systems in
[1]. To obtain canonical conserved densities, a diagonalisation procedure was
used. This allows to split the matrix relation (15) into $m$ scalar relations
similar to (7) and to apply formula (8). Unfortunately, many of the known
integrable evolution systems do not satisfy the non-degeneracy condition. In
particular, systems (2) are not non-degenerate.

In this paper we generalise the symmetry approach to a class of systems (3). A
fundamental point in our work is that for systems (3), albeit being degenerate
and non-polynomial, it is possible to perform a full diagonalisation of the
Fréchet derivative in the jet space. As a result, a formal recursion operator
also reduces to the diagonal form. This allows us to generate explicit
integrability conditions in the form of canonical conservation laws. Using
these conditions we perform a partial classification of systems of the form

(18) $\displaystyle u_{t}$ $\displaystyle=v,$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+f(u,u_{1},\ldots,u_{3},v,v_{1})$

without the polinomiality restriction for $f$. Most of the obtained equations
are not polynomial.

This paper is organised as follows. Section 2 is devoted to the
diagonalisation procedure of systems (3). In Subsection 2.1 we recall results
from [1] related to the diagonalisation procedure for non-degenerate systems.
In Subsection 2.2 we introduce the more general concept of regularly
diagonalisable systems. We formulate all statements about such systems without
proofs, because they can be proved exactly in the same way as the
corresponding statements for the non-degenerate case (cf. [1]).

Section 3 is devoted to the integrability conditions for systems (3). In
Section 4 we present, in the case $k=2$, explicit formulas for several first
coefficients of relevant matrix pseudo-differential series appearing in the
diagonalisation, as well as the simplest canonical densities.

In Section 5 we show that any $S$-integrable system (18) has the form

(19) $\displaystyle u_{t}$ $\displaystyle=v,$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+g_{u_{2}}u_{3}^{2}+2g_{v}v_{1}u_{3}+2\left(g_{u_{1}}u_{2}+g_{u}u_{1}\right)u_{3}+f_{2}v_{1}^{2}+f_{1}v_{1}+f_{0},$

where $g=g(u,u_{1},u_{2},v)$ and $f_{i}=f_{i}(u,u_{1},u_{2},v)$.

In Section 6 we find all $S$-integrable systems (19) with $g=g(u)$ while
Section 7 is devoted to classification of integrable systems with
$g=g(u_{1})$. In these sections some finite lists of systems are found.

Since canonical densities provide only necessary integrability conditions,
actual integrability for obtained systems should be justified independently.
One of the possible ways to do that is to construct a compact form for the
formal recursion operators, whose first coefficients we initially found in
Sections 5 and 6. The main tool here is using quasi-local anzatzes for
recursion, Hamiltonian and symplectic operators (see [14, 15, 16, 17]).
Written in this form, a formal recursion operator is a usual recursion
operator, and it can be applied to some simple “seed” symmetries to generate
the hierarchy of infinitesimal symmetries. In Section 8 we demonstrate this
technique applied to the systems found in Section 5.

##  2\. diagonalisation

###  2.1. Non-degenerate systems

For non-degenerate systems (13) (see Definition 2), integrability conditions
can be obtained through the diagonalisation of the operator
$D_{t}-\mathbf{\Phi}_{*}$, which corresponds to the linearisation of system
(13) around an arbitrary solution ${\bf u}$. Namely, one can show [1,
Proposition 2.1] that there exists a matrix pseudo-differential series

(20) $\mathbf{T}=\mathbf{T}_{0}+\sum_{l>0}\mathbf{T}_{-l}D^{-l},$

where $\mathbf{T}_{i}$ are $m\times m$ matrices depending on $\mathbf{u}$ and
its derivatives, such that

(21)
$\mathbf{T}^{-1}(D_{t}-\mathbf{\Phi}_{*})\mathbf{T}=D_{t}-\overline{\mathbf{\Phi}}_{*},$

where
$\overline{\mathbf{\Phi}}_{*}=\mathbf{T}^{-1}\mathbf{\Phi}_{*}\mathbf{T}-\mathbf{T}^{-1}\mathbf{T}_{t}$
is a diagonal matrix pseudo-differential series

(22) $\overline{\mathbf{\Phi}}_{*}={\rm
diag}(\Phi_{1},\dots,\Phi_{m}),\qquad\Phi_{i}=\sum_{l=-n}^{\infty}p_{-l,i}D^{-l},$

such that $p_{n,i}\neq p_{n,j}$ for $i\neq j$ and $p_{n,i}\neq 0$ for all
possible $i$. In (20) $\mathbf{T}_{0}$ is a matrix that diagonalises the
separant matrix $\mathbf{\Sigma}$ of system (13):

(23)
$\boldsymbol{\Lambda}=\mathbf{T}_{0}^{-1}\,\mathbf{\Sigma}\,\mathbf{T}_{0}={\rm
diag}(p_{n,1},\dots,p_{n,m}),$

and

(24)
$\overline{\mathbf{\Phi}}_{*}=\boldsymbol{\Lambda}D^{n}+\sum_{l=-n+1}^{\infty}\overline{\mathbf{\Phi}}^{(-l)}D^{-l}.$

After diagonalisation, formal recursion and symplectic operators satisfy the
equations

(25)
$\displaystyle\overline{\mathbf{R}}_{t}=[\overline{\mathbf{\Phi}}_{*},\overline{\mathbf{R}}],$
(26)
$\displaystyle\overline{\mathbf{S}}_{t}+\overline{\mathbf{S}}\,\overline{\mathbf{\Phi}}_{*}+\overline{\mathbf{\Phi}}_{*}^{+}\,\overline{\mathbf{S}}=0,$

where

$\overline{\mathbf{R}}=\mathbf{T}^{-1}\mathbf{R}\mathbf{T},\qquad\overline{\mathbf{S}}=\mathbf{T}^{+}\mathbf{S}\mathbf{T}.$

The following result can be easily deduced (see [1, Theorem 2.1]) from (25):

######  Proposition 1.

If a non-degenerate system (13) possesses a formal recursion operator
$\mathbf{R}$, then the operator $\overline{\mathbf{R}}$ is diagonal:

(27) $\overline{\mathbf{R}}={\rm diag}(R_{1},\dots,R_{m}),\qquad
R_{i}=\sum_{l=s_{i}}^{\infty}r_{-l,i}D^{-l},\qquad r_{-s_{i},i}\neq 0,\
s_{i}\in\mathbb{Z}.$

As both $\overline{\mathbf{\Phi}}_{*}$ and $\overline{\mathbf{R}}$ are
diagonal, each of the series $R_{i}$ satisfies a scalar relation

(28) $(R_{i})_{t}=[\Phi_{i},\,R_{i}].$

We will call a formal recursion operator $\mathbf{R}$ _non-degenerate_ if in
formula (27) we have $s_{i}\neq 0$ for all $i$. Without loss of generality one
may assume that in the non-degenerate case the numbers $s_{i}$ are equal to
$-1$. Just as in the scalar case, the functions
$\rho_{ij}={\operatorname{res}}\,R_{i}^{j}$ are local conserved densities for
(13).

The existence a non-degenerate formal recursion operator can be chosen as an
integrability criterion for a non-degenerate evolution system (13). Without
imposing the nondegeneracy condition to $\mathbf{R}$, the following problem
arises. Suppose that the system consists of two uncoupled systems, one
integrable and the other non-integrable. Then the whole system is not
integrable but nevertheless it has a degenerate formal recursion operator.

Similar problems arise on the level of symmetries or conservation laws: for
multi-component systems it is not enough for integrability to have one
infinite sequence of symmetries or conservation laws. These sequences must be
non-degenerate in some sense.

A symmetry

(29) ${\bf u}_{\tau}=\boldsymbol{\psi}({\bf u},{\bf u}_{1},\dots,{\bf u}_{k})$

is called _non-singular_ if the separant matrix of $\boldsymbol{\psi}$ is non-
singular [1, p. 10]. It can be shown that if the system has an infinite
sequence of non-singular symmetries then it possesses a non-degenerate formal
recursion operator.

The requirement for nonsingularity of symmetries can be weakened. For
instance, the following statement can be formulated using the diagonalisation
procedure. It is easy to see that the conjugation by $\mathbf{T}$ diagonalises
not only $\mathbf{\Phi}_{*}$ and $\mathbf{R}$, but also the Fréchet derivative
$\mathbf{\Psi}_{*}$ of any symmetry of the system (13). Namely, we have

$\mathbf{T}^{-1}(D_{\tau}-\mathbf{\Psi}_{*})\mathbf{T}=D_{\tau}-\overline{\mathbf{\Psi}}_{*},$

where $\overline{\mathbf{\Psi}}_{*}$ is diagonal:

$\overline{\mathbf{\Psi}}_{*}={\rm
diag}(\Psi_{1},\dots,\Psi_{m}),\qquad\Psi_{i}=\sum_{l=-d_{i}}^{\infty}p_{-l,i}D^{-l},\qquad
p_{d_{i},i}\neq 0.$

The sequence of symmetries

$\mathbf{u}_{\tau_{j}}=\boldsymbol{\psi}_{j},\qquad j\to\infty$

is called non-degenerate if $(d_{k})_{j}\to\infty$ for any $1\leq k\leq m$.

The following statement can be proved in the same way as Theorem 1.7 from
[11]:

######  Proposition 2.

Suppose that a non-degenerate system (13) has a non-degenerate sequence of
symmetries. Then the system has a non-degenerate formal recursion operator.

According to Proposition 1, any formal recursion operator becomes diagonal
when the Fréchet derivative of the system is already diagonalised.

The form of symplectic operators of non-degenerate systems (13) depends on the
order $n$ of the series $\mathbf{\Phi}_{*}$ and the number of equations $m$ in
system (13) (see [1]). Write the series $\overline{\mathbf{S}}$ in (26) as
$\overline{\mathbf{S}}=\sum_{i=-l}^{\infty}\overline{\mathbf{S}}^{(-i)}D^{-i}$.
The coefficient of the highest power of $D$, $D^{l+n}$, in (26) implies that

(30)
$(-1)^{n}\boldsymbol{\Lambda}\overline{\mathbf{S}}^{(l)}+\overline{\mathbf{S}}^{(l)}\boldsymbol{\Lambda}=0,$

where $\boldsymbol{\Lambda}$ is defined by (23). The equation (26) can be
expanded as

$\sum_{k\geq
0}\left[(-1)^{n}\boldsymbol{\Lambda}\overline{\mathbf{S}}^{(l-k)}+\overline{\mathbf{S}}^{(l-k)}\boldsymbol{\Lambda}\right]D^{l+n-k}+\sum_{i\geq
0}\overline{\mathbf{S}}^{(l-i)}_{t}D^{l-i}\\\ +\sum_{k\geq
1}\sum_{i=0}^{k-1}\sum_{j=0}^{k-i}\Bigg{[}(-1)^{n-j}\binom{n-j}{k-j-i}D^{k-i-j}\left(\overline{\mathbf{\Phi}}^{(n-j)}\overline{\mathbf{S}}^{(l-i)}\right)\\\
{}+\overline{\mathbf{S}}^{(l-i)}\binom{l-i}{k-j-i}D^{k-j-i}\left(\overline{\mathbf{\Phi}}^{(n-j)}\right)\Bigg{]}D^{l+n-k}=0.$

We see that equations from lower coefficients in (26) are of the form

(31)
$(-1)^{n}\boldsymbol{\Lambda}\overline{\mathbf{S}}^{(l-k)}+\overline{\mathbf{S}}^{(l-k)}\boldsymbol{\Lambda}+\boldsymbol{\Delta}_{k}=0$

where $\boldsymbol{\Delta}_{k}$ depends on previously found
$\overline{\mathbf{S}}^{(i)}$’s and on $\overline{\mathbf{\Phi}}^{(i)}$. The
terms in $\boldsymbol{\Delta}_{k}$ are products
$D^{i}\overline{\mathbf{\Phi}}^{(j)}\cdot D^{r}\overline{\mathbf{S}}^{(s)}$ or
$\overline{\mathbf{S}}^{(i)}$, with factors $\overline{\mathbf{S}}^{(i)}$
previously found.

  1. (1)

If $n$ is odd, (30) implies that $\overline{\mathbf{S}}^{(l)}$ is diagonal.
Inductively, $\boldsymbol{\Delta}_{k}$ is also diagonal, which implies that,
$\overline{\mathbf{S}}^{(l-k)}$ must be diagonal, and an integrability
condition $\boldsymbol{\Delta}_{k}=0$. Thus we can write

(32) $\overline{\mathbf{S}}={\rm diag}(S_{1},\dots,S_{m}),\qquad
S_{i}=\sum_{j=-k_{i}}^{\infty}s_{-j,i}D^{-j}.$

We will call a formal symplectic operator $\overline{\mathbf{S}}$ in formula
(32) _non-degenerate_ if $k_{i}\neq 0$ for all $i$. For systems that have both
non-degenerate formal recursion and symplectic operators we may assume, as in
the scalar case, that $k_{i}=1$ and that

(33)
$\overline{\mathbf{R}}^{+}=-\overline{\mathbf{S}}\,\overline{\mathbf{R}}\,\overline{\mathbf{S}}^{-1},\qquad\overline{\mathbf{S}}^{+}=-\overline{\mathbf{S}}.$

  2. (2)

If $n$ is even, $n=2k$, it follows from (30) that
$(\overline{\mathbf{S}}_{l})^{i}_{j}\neq 0$ only if $p_{n,i}+p_{n,j}=0$.
Rearrange the entries of $\overline{\mathbf{\Psi}}_{*}$ to obtain a matrix
$\boldsymbol{\Lambda}$ of the form

(34)
$\operatorname{diag}(p_{n,1},-p_{n,1},\cdots,p_{n,s},-p_{n,s},p_{n,2s+1},\cdots,p_{n,2k}),$

where $p_{n,i}+p_{n,j}\neq 0$ for $i,j>2s$. Formulas (31) then imply that
$\overline{\mathbf{S}}$ has the $2\times 2$ block form

$\overline{\mathbf{S}}={\rm
diag}(\mathbf{S}_{1},\dots,\mathbf{S}_{s},\mathbf{0},\dots,\mathbf{0}),\qquad\mathbf{S}_{i}=\begin{bmatrix}0&\sigma_{i}\\\
\tau_{i}&0\end{bmatrix}$

with $\sigma_{i},\tau_{i}$ being scalar pseudodifferential series. If the
system has a non-degenerate formal recursion operator, then each non-constant
block $\mathbf{S}_{i}$ can be reduced (see Remark 2) to

(35) $\mathbf{S}_{i}=\begin{bmatrix}0&\sigma_{i}\\\
-\sigma_{i}^{+}&0\end{bmatrix},\qquad\sigma_{i}=\sum_{j=-1}^{\infty}s_{-j,i}D^{-j}.$

A formal symplectic operator is called non-degenerate if it has the form

$\overline{\mathbf{S}}={\rm diag}(\mathbf{S}_{1},\dots,\mathbf{S}_{k}),$

where the blocks $\mathbf{S}_{i}$ are of the form (35). In this case, if $p$
is an eigenvalue of $\boldsymbol{\Lambda}$, then $-p$ is also an eigenvalue.

######  Definition 3.

Non-degenerate systems (13) that have both non-degenerate formal recursion and
symplectic operators are called S-integrable.

######  Remark 3.

The existence of an infinite sequence of non-singular222This means that the
leading coefficient of the differential operator
$(\delta\rho/\delta\mathbf{u})_{*}$ is non-singular [1, p. 10]. conservation
laws

$D_{t}\mathbf{\rho}=D\sigma$

implies (cf. [11]) the existence of both non-degenerate formal recursion and
symplectic operators.

###  2.2. Regularly diagonalisable systems

A similar diagonalisation is possible for some degenerate systems (13).
Namely, it can happen that there exists a pseudo-differential matrix series
$\mathbf{T}$ such that the series $\overline{\mathbf{\Phi}}_{*}$ in (21) takes
the form

$\overline{\mathbf{\Phi}}_{*}={\rm
diag}(\Phi_{1},\dots,\Phi_{m}),\qquad\Phi_{i}=\sum_{l=-k}^{\infty}p_{-l,i}D^{-l},$

where $1<k<n$, $p_{k,i}\neq p_{k,j}$ for $i\neq j$ and $p_{k,i}\neq 0$ for all
$i$333Most statements can be easily generalised to the case when one of the
$p_{k,i}$ is equal to zero.. In this case we say that (13) admits a  _regular
diagonalisation of order_ $k$. In [1, Proposition 2.1] it was shown that non-
degenerate systems (13) are regularly diagonalisable of order $n$.

Our main observation is that systems (3) are regularly diagonalisable.

All the statements from [1] about non-degenerate systems reviewed in the
previous Section, can be easily generalised to the case of systems admitting a
regular diagonalisation. For such systems, Fréchet derivatives of symmetries
and formal recursion operators become diagonal simultaneously with the Fréchet
derivative of the system. The definition of the order of a symmetry is the
same to that of non-degenerate systems. The canonical form of a formal
symplectic operator depends on whether the number $k$ is even or odd.

Recursion and symplectic operators. Following the line in [1, Theorem 2.1] it
can be proved that, for regularly diagonalisable $S$-integrable systems (see
Definition 1) and without loss of generality, formal recursion and symplectic
operators $\overline{\mathbf{R}}$ and $\overline{\mathbf{S}}$ can be assumed
to have order $1$ and be related by the relations (33)

$\overline{\mathbf{R}}^{+}=-\overline{\mathbf{S}}\,\overline{\mathbf{R}}\,\overline{\mathbf{S}}^{-1},\qquad\overline{\mathbf{S}}^{+}=-\overline{\mathbf{S}}.$

The matrix series $\overline{\mathbf{R}}$ is always diagonal for regularly
diagonalisable systems. The series $\overline{\mathbf{S}}$ is diagonal if $k$
is odd and block-diagonal if $k$ is even, similarly to the non-degenerate case
described in Section 2.1.

Integrability conditions. After diagonalisation, relation (15) splits into $m$
scalar equations (28). These equations define scalar canonical conservation
laws

(36) $(\rho_{i,j})_{t}=(\sigma_{i,j})_{x},\qquad i=1,\dots,m,\quad
j=-1,0,1,2,\dots,$

where

$\rho_{i,j}={\rm res}\,R_{i}^{j}\quad{\rm for}\,\,j\neq 0\qquad{\rm
and}\quad\rho_{i,0}={\rm res}\,{\rm log}(R_{i}).$

As in the scalar case (see Remark 2), some linear combinations of canonical
conserved densities have to be trivial (i.e. they are total $x$-derivatives).
In the next section we describe in detail these integrability conditions in
the case of systems (3).

##  3\. Systems of the Boussinesq type

Consider systems of the form (3). Their Fréchet derivative and separant are

(37) ${\mathbf{\Phi}}_{*}=\begin{bmatrix}[c]0&1\\\
U&V\end{bmatrix},\qquad\mathbf{\Sigma}=\begin{bmatrix}[c]0&0\\\
f_{u_{2k}}&0\end{bmatrix},$

where, denoting $f_{u_{i}}=\frac{\partial f}{\partial u_{i}},$
$f_{v_{i}}=\frac{\partial f}{\partial v_{i}}$

$U\stackrel{{\scriptstyle\text{def}}}{{=}}f_{u_{2k}}D^{2k}+\cdots+f_{u},\qquad
V\stackrel{{\scriptstyle\text{def}}}{{=}}f_{v_{k-1}}D^{k-1}+\cdots+f_{v}.$

######  Theorem 1.

For any system (3) with $k>1$ there exists a matrix pseudo-differential series
$\mathbf{T}$ of the form

(38) $\mathbf{T}=\begin{bmatrix}[c]1&T_{2}^{-1}\\\ T_{1}&1\end{bmatrix},$

being $T_{1}$ and $T_{2}$ scalar pseudo-differential series of order $k$, such
that

$\mathbf{T}^{-1}(D_{t}-\mathbf{\Phi}_{*})\mathbf{T}=D_{t}-\overline{\mathbf{\Phi}}_{*},$

where

(39)
$\overline{\mathbf{\Phi}}_{*}=\begin{bmatrix}f_{u_{2k}}^{1/2}D^{k}+\cdots&0\\\
0&-f_{u_{2k}}^{1/2}D^{k}+\cdots\end{bmatrix}.$

###### Proof.

Substituting

$\overline{\mathbf{\Phi}}_{*}=\begin{bmatrix}[c]\Phi_{1}&0\\\
0&\Phi_{2}\end{bmatrix},\qquad\mathbf{T}=\begin{bmatrix}[c]1&T_{2}^{-1}\\\
T_{1}&1\end{bmatrix}$

into

(40)
$0=\mathbf{T}\overline{\mathbf{\Phi}}_{*}-\mathbf{\Phi}_{*}\mathbf{T}+\mathbf{T}_{t},$

we obtain that both $T_{1}$ and $T_{2}$ satisfy the Riccati-type operator
equation

(41) $T_{t}+T^{2}-VT=U.$

Let us describe the solutions to this equation writing them in the form

$T=\tau_{k}D^{k}+\tau_{k-1}D^{k-1}+\cdots+\tau_{0}+\tau_{-1}D^{-1}+\cdots.$

The coefficient of $D^{2k}$ in (41) implies that $\tau_{k}^{2}=f_{u_{2k}}$.
Relations which follow from lower powers of $D$ have the form

$2\tau_{k}\tau_{j}=c(\tau_{k},\tau_{k-1},\ldots,\tau_{j+1}),\quad j<k$

where the rhs $c$ is a differential expression depending only on previously
found unknowns $\tau_{i}$, $i>j$ and the coefficients of $U$ and $V$, i.e. on
the rhs of system (3). Thus, choosing a leading coefficient $\tau_{k}$, we
uniquely find all other coefficients of $T$. As a result, we get two different
solutions of equation (41) with leading coefficients $f_{u_{2k}}^{1/2}$ and
$-f_{u_{2k}}^{1/2}$. Taking the first solution as $T_{1}$ and the second as
$T_{2}$ leads to formula (39). ∎

######  Remark 4.

The diagonalisation procedure of Theorem 1 can be applied to many more
systems. For example, the systems in Remark 1 can be diagonalised using a
matrix (38) with $T_{1}$ and $T_{2}$ being two different solutions of the
operator equation

$T_{t}+TD\cdot T-VT=U.$

The resulting diagonalised Fréchet derivative takes the form

$\overline{\mathbf{\Phi}}_{*}=\begin{bmatrix}f_{u_{2k-1}}^{1/2}D^{k}+\cdots&0\\\
0&-f_{u_{2k-1}}^{1/2}D^{k}+\cdots\end{bmatrix}$

and, thus, this class of systems is regularly diagonalisable.

Recursion, symplectic operators, and integrability conditions. As it was
mentioned in Paragraphs 2.1 and 2.2, in the case $k=2l$ the matrix
$\overline{\mathbf{S}}$ consists of $l$ blocks of size $2\times 2$. In
particular, for systems (3) with $k=2$ we have

(42) $\overline{\mathbf{\Phi}}_{*}=\begin{bmatrix}[c]\Phi_{1}&0\\\
0&\Phi_{2}\end{bmatrix},\qquad\overline{\mathbf{R}}=\begin{bmatrix}[c]R_{1}&0\\\
0&R_{2}\end{bmatrix},\qquad\overline{\mathbf{S}}=\begin{bmatrix}[c]0&S_{1}\\\
S_{2}&0\end{bmatrix}.$

Considering (33) it follows that

(43) $\overline{\mathbf{R}}=\begin{bmatrix}[c]R&0\\\
0&-S^{-1}R^{+}S\end{bmatrix},\qquad\overline{\mathbf{S}}=\begin{bmatrix}[c]0&S\\\
-S^{+}&0\end{bmatrix}.$

where relations (25), (26) imply that

(44) $R_{t}=[\Phi_{1},R],\qquad S_{t}+\Phi_{1}^{+}S+S\Phi_{2}=0.$

Using (43) it is easy to prove that

(45) $\rho_{1,j}+(-1)^{j}\rho_{2,j}\in{\rm Im}\,D,\qquad j=-1,0,1,\dots.$

Let us denote

$\displaystyle\rho_{-1}=\tfrac{1}{2}(\rho_{1,-1}+\rho_{2,-1}),$
$\displaystyle\begin{gathered}\quad\rho_{4i}=\rho_{1,2i}+\rho_{2,2i},\quad\rho_{4i+1}=\rho_{1,2i}-\rho_{2,2i},\\\
\rho_{4i+2}=\rho_{1,2i+1}-\rho_{2,2i+1},\quad\rho_{4i+3}=\rho_{1,2i+1}+\rho_{2,2i+1},\end{gathered}\qquad
i=0,1,2,3,\ldots.$

According to (36), all functions $\rho_{i}$ have to be local conserved
densities: $D_{t}\rho_{i}=D\sigma_{i}$ for $i=-1,0,1,\ldots$. From (45) it
follows that the even-numbered densities $\rho_{2i}$ are divergencies (total
$x$-derivatives): $\rho_{2i}=D\pi_{2i}$, $i=0,1,2,\ldots$. In this case
$\sigma_{2i}=D_{t}\pi_{2i}$. The fluxes $\sigma_{2i-1}$ and potentials
$\pi_{2i}$ are functions appearing in the coefficients of $R$ and $S$.

The relations

(46) $\displaystyle D_{t}\rho_{i}$ $\displaystyle=D\sigma_{i},$ if $i$ is $-1$
or odd $\displaystyle\rho_{i}$ $\displaystyle=D\pi_{i},$ if $i$ is even

form an infinite set of integrability conditions for $S$-integrable systems
(3), $k=2$.

##  4\. Boussinesq-type systems of fourth order

Let us detail the results and constructions of the previous section in the
case of systems (3), $k=2$:

(47) $u_{t}=v,\qquad v_{t}=f(u,u_{1},u_{2},u_{3},u_{4},v,v_{1}).$

For such systems $U$ and $V$ are differential operators of orders 4 and 1:

$U=\mathfrak{u}_{4}D^{4}+\mathfrak{u}_{3}D^{3}+\mathfrak{u}_{2}D^{2}+\mathfrak{u}_{1}D+\mathfrak{u}_{0},\qquad
V=\mathfrak{v}_{1}D+\mathfrak{v}_{0}$

where for conciseness we have denoted $\mathfrak{u}_{i}=\partial f/\partial
u_{i}$, $\mathfrak{u}_{0}=\partial f/\partial u$, $\mathfrak{v}_{1}=\partial
f/\partial v_{1}$ and $\mathfrak{v}_{0}=\partial f/\partial v$.

Solving (41) yields that the first terms in $T_{1}$ are

$T_{1}=\sqrt{\mathfrak{u}_{4}}\,D^{2}+\left[\frac{\mathfrak{u}_{3}}{2\sqrt{\mathfrak{u}_{4}}}+\frac{\mathfrak{v}_{1}}{2}-\frac{D_{x}\mathfrak{u}_{4}}{2\sqrt{\mathfrak{u}_{4}}}\right]\,D\\\
+\bigg{[}-\frac{\mathfrak{u}_{3}^{2}}{8\mathfrak{u}_{4}^{3/2}}+\frac{3D\mathfrak{u}_{4}\,\mathfrak{u}_{3}}{8\mathfrak{u}_{4}^{3/2}}+\frac{\mathfrak{v}_{1}^{2}}{8\sqrt{\mathfrak{u}_{4}}}-\frac{\left(D\mathfrak{u}_{4}\right)^{2}}{8\mathfrak{u}_{4}^{3/2}}+\frac{\mathfrak{v}_{0}}{2}-\frac{D\mathfrak{u}_{3}}{2\sqrt{\mathfrak{u}_{4}}}\\\
\quad+\frac{\mathfrak{v}_{1}D\mathfrak{u}_{4}}{8\mathfrak{u}_{4}}-\frac{1}{2}D\mathfrak{v}_{1}-\frac{D_{t}\mathfrak{u}_{4}}{4\mathfrak{u}_{4}}+\frac{D^{2}\mathfrak{u}_{4}}{4\sqrt{\mathfrak{u}_{4}}}+\frac{\mathfrak{u}_{2}}{2\sqrt{\mathfrak{u}_{4}}}\bigg{]}{}+\cdots$

and $T_{2}$ can be obtained from $T_{1}$ via the replacement
$\sqrt{\mathfrak{u}_{4}}\to-\sqrt{\mathfrak{u}_{4}}$. Relation (40) leads to

$\Phi_{1}=\sqrt{\mathfrak{u}_{4}}\,D^{2}+\bigg{[}\frac{\mathfrak{u}_{3}}{2\sqrt{\mathfrak{u}_{4}}}+\frac{\mathfrak{v}_{1}}{2}-\frac{D\mathfrak{u}_{4}}{2\sqrt{\mathfrak{u}_{4}}}\bigg{]}\,D\\\
{}+\bigg{[}-\frac{\mathfrak{u}_{3}^{2}}{8\mathfrak{u}_{4}^{3/2}}+\frac{3D\mathfrak{u}_{4}\mathfrak{u}_{3}}{8\mathfrak{u}_{4}^{3/2}}+\frac{\mathfrak{v}_{1}^{2}}{8\sqrt{\mathfrak{u}_{4}}}-\frac{\left(D\mathfrak{u}_{4}\right)^{2}}{8\mathfrak{u}_{4}^{3/2}}+\frac{\mathfrak{v}_{0}}{2}-\frac{D\mathfrak{u}_{3}}{2\sqrt{\mathfrak{u}_{4}}}\\\
\quad{}+\frac{\mathfrak{v}_{1}D\mathfrak{u}_{4}}{8\mathfrak{u}_{4}}-\frac{1}{2}D\mathfrak{v}_{1}-\frac{D_{t}\mathfrak{u}_{4}}{4\mathfrak{u}_{4}}+\frac{D_{xx}\mathfrak{u}_{4}}{4\sqrt{\mathfrak{u}_{4}}}+\frac{\mathfrak{u}_{2}}{2\sqrt{\mathfrak{u}_{4}}}\bigg{]}+\cdots$

The series $\Phi_{2}$ is related to $\Phi_{1}$ by the replacement
$\sqrt{\mathfrak{u}_{4}}\to-\sqrt{\mathfrak{u}_{4}}$.

The first terms of series $R$ and $S$ defined by relations (44) have the form

$R=\mathfrak{u}_{4}^{1/4}\,D+\left[-\frac{1}{2}\sigma_{-1}+\frac{\mathfrak{v}_{1}}{4\mathfrak{u}_{4}^{1/4}}-\frac{3D\mathfrak{u}_{4}}{8\mathfrak{u}_{4}^{3/4}}+\frac{\mathfrak{u}_{3}}{4\mathfrak{u}_{4}^{3/4}}\right]\\\
{}+\Bigg{[}\frac{\sigma_{-1}^{2}}{8\mathfrak{u}_{4}^{1/4}}+\frac{1}{4}D\sigma_{-1}+\frac{\sigma_{0}}{2\mathfrak{u}_{4}^{1/4}}+\frac{\mathfrak{v}_{1}^{2}}{32\mathfrak{u}_{4}^{3/4}}-\frac{35\left(D\mathfrak{u}_{4}\right)^{2}}{128\mathfrak{u}_{4}^{7/4}}+\frac{\mathfrak{v}_{0}}{4\mathfrak{u}_{4}^{1/4}}-\frac{\mathfrak{u}_{3}\mathfrak{v}_{1}}{16\mathfrak{u}_{4}^{5/4}}-\frac{3D\mathfrak{u}_{3}}{8\mathfrak{u}_{4}^{3/4}}\\\
{}+\frac{3\mathfrak{v}_{1}D\mathfrak{u}_{4}}{16\mathfrak{u}_{4}^{5/4}}+\frac{3\mathfrak{u}_{3}D\mathfrak{u}_{4}}{8\mathfrak{u}_{4}^{7/4}}-\frac{3D\mathfrak{v}_{1}}{8\mathfrak{u}_{4}^{1/4}}-\frac{D_{t}\mathfrak{u}_{4}}{8\mathfrak{u}_{4}^{5/4}}+\frac{5D^{2}\mathfrak{u}_{4}}{16\mathfrak{u}_{4}^{3/4}}+\frac{\mathfrak{u}_{2}}{4\mathfrak{u}_{4}^{3/4}}-\frac{3\mathfrak{u}_{3}^{2}}{32\mathfrak{u}_{4}^{7/4}}\Bigg{]}\,D^{-1}+\cdots$

and

(48) $S={\rm e}^{\pi_{0}}D+{\rm
e}^{\pi_{0}}\left[-\frac{\pi_{2}}{\mathfrak{u}_{4}^{1/4}}-\frac{\mathfrak{v}_{1}}{4\sqrt{\mathfrak{u}_{4}}}-\frac{3D\mathfrak{u}_{4}}{8\mathfrak{u}_{4}}+\frac{\sigma_{-1}}{2\mathfrak{u}_{4}^{1/4}}+\frac{\mathfrak{u}_{3}}{4\mathfrak{u}_{4}}\right]\\\
{}+{\rm
e}^{\pi_{0}}\Bigg{[}\frac{\pi_{4}}{2\sqrt{\mathfrak{u}_{4}}}+\frac{\pi_{2}^{2}}{2\sqrt{\mathfrak{u}_{4}}}+\frac{\sigma_{1}}{4\sqrt{\mathfrak{u}_{4}}}+\frac{\sigma_{-1}^{2}}{8\sqrt{\mathfrak{u}_{4}}}+\frac{\mathfrak{v}_{1}^{2}}{32\mathfrak{u}_{4}}-\frac{35(D\mathfrak{u}_{4})^{2}}{128\mathfrak{u}_{4}^{2}}-\frac{3D\mathfrak{u}_{3}}{8\mathfrak{u}_{4}}+\frac{3\mathfrak{u}_{3}D\mathfrak{u}_{4}}{8\mathfrak{u}_{4}^{2}}\\\
{}-\frac{D_{t}\mathfrak{u}_{4}}{8\mathfrak{u}_{4}^{3/2}}+\frac{5D^{2}\mathfrak{u}_{4}}{16\mathfrak{u}_{4}}+\frac{\mathfrak{u}_{2}}{4\mathfrak{u}_{4}}-\frac{3\mathfrak{u}_{3}^{2}}{32\mathfrak{u}_{4}^{2}}\Bigg{]}D^{-1}+\cdots$

where the functions $\sigma_{2i-1}$ and $\pi_{2i}$, $i=0,1,2,\ldots$ are those
in (46).

The first canonical densities444The expressions for the densities $\rho_{i}$
are defined up to a total derivative $\rho_{i}\to\rho_{i}+D\alpha_{i}$; here
we choose $\alpha_{i}$ to get simpler expressions for the densities, without
forgetting to adjust the fluxes correspondingly,
$\sigma_{i}\to\sigma_{i}+D_{t}\alpha_{i}$. take the form

(49) $\displaystyle\rho_{-1}=\frac{1}{\mathfrak{u}^{1/4}},$ (50)
$\displaystyle\rho_{0}=\frac{1}{2}\frac{\mathfrak{u}_{3}}{\mathfrak{u}_{4}}-\frac{3D\mathfrak{u}_{4}}{4\mathfrak{u}_{4}},$
(51)
$\displaystyle\rho_{1}=-\frac{\sigma_{-1}}{4\mathfrak{u}_{4}}+\frac{1}{2}\frac{\mathfrak{v}_{1}}{\sqrt{\mathfrak{u}_{4}}},$
(52)
$\displaystyle\rho_{2}=\frac{\sigma_{0}}{2\mathfrak{u}_{4}^{1/4}}+\frac{1}{2}D\sigma_{-1}+\frac{3\mathfrak{v}_{1}D\mathfrak{u}_{4}}{8\mathfrak{u}_{4}^{5/4}}-\frac{3D\mathfrak{v}_{1}}{4\mathfrak{u}_{4}^{1/4}}+\frac{\mathfrak{v}_{0}}{2\mathfrak{u}_{4}^{1/4}}-\frac{\mathfrak{u}_{3}\mathfrak{v}_{1}}{8\mathfrak{u}_{4}^{5/4}},$
(53) $\displaystyle\begin{aligned}
\rho_{3}=&\frac{\sigma_{1}}{2\mathfrak{u}_{4}^{1/4}}+\frac{\sigma_{-1}^{2}}{4\mathfrak{u}_{4}^{1/4}}-\frac{D_{t}\mathfrak{u}_{4}}{4\mathfrak{u}_{4}^{5/4}}-\frac{35\left(D\mathfrak{u}_{4}\right)^{2}}{64\mathfrak{u}_{4}^{7/4}}-\frac{3D\mathfrak{u}_{3}}{4\mathfrak{u}_{4}^{3/4}}\\\
&{}+\frac{3\mathfrak{u}_{3}D\mathfrak{u}_{4}}{4\mathfrak{u}_{4}^{7/4}}+\frac{5D^{2}\mathfrak{u}_{4}}{8\mathfrak{u}_{4}^{3/4}}-\frac{3\mathfrak{u}_{3}^{2}}{16\mathfrak{u}_{4}^{7/4}}+\frac{\mathfrak{v}_{1}^{2}}{16\mathfrak{u}_{4}^{3/4}}+\frac{\mathfrak{u}_{2}}{2\mathfrak{u}_{4}^{3/4}}.\end{aligned}$

Subsequent densities can be calculated by computer. Their length grows fast,
and we do not present them here.

##  5\. On the classification of integrable fourth-order Boussinesq-type
equations with separant 1

As an application of the previous developments, we consider systems (47) of
the form

(54) $u_{t}=v,\qquad v_{t}=u_{4}+f(u,u_{1},u_{2},u_{3},v,v_{1})$

and show how to use the integrability conditions for finding the integrable
cases.

For systems (54), $\mathfrak{u}_{4}$ is equal to 1, considerably simplifying
the canonical densities, being the first ones:

$\displaystyle\rho_{-1}=1,\quad\rho_{0}=\frac{\mathfrak{u}_{3}}{2},\quad\rho_{1}=\frac{\mathfrak{v}_{1}}{2},\quad\rho_{2}=\frac{\sigma_{0}}{2}-\frac{3}{4}D\mathfrak{v}_{1}-\frac{\mathfrak{u}_{3}\mathfrak{v}_{1}}{8}+\frac{\mathfrak{v}_{0}}{2},$
$\displaystyle\rho_{3}=\frac{\sigma_{1}}{2}-\frac{3}{4}D\mathfrak{u}_{3}-\frac{3\mathfrak{u}_{3}^{2}}{16}+\frac{\mathfrak{v}_{1}^{2}}{16}+\frac{\mathfrak{u}_{2}}{2},$
$\displaystyle\begin{split}\rho_{4}&=\sigma_{2}-\frac{1}{2}D_{t}\mathfrak{u}_{3}+\frac{3}{4}\mathfrak{u}_{3}D\mathfrak{u}_{3}-D\mathfrak{u}_{2}+\frac{1}{2}D^{2}\mathfrak{u}_{3}-\frac{1}{4}\mathfrak{v}_{1}D\mathfrak{v}_{1}\\\
&\qquad{}-\frac{1}{8}\mathfrak{u}_{3}\mathfrak{v}_{1}^{2}+\frac{\mathfrak{u}_{3}^{3}}{8}+\frac{\mathfrak{v}_{0}\mathfrak{v}_{1}}{2}-\frac{\mathfrak{u}_{2}\mathfrak{u}_{3}}{2}+\mathfrak{u}_{1}.\end{split}$

We will refer to the integrability condition corresponding to a density
$\rho_{i}$ as $C_{i}$.

######  Lemma 1.

If system (54) satisfies the first two integrability conditions

(55) $\displaystyle C_{0}:\,\mathfrak{u}_{3}=2Dg(u,u_{1},u_{2},v)$ (56)
$\displaystyle C_{1}:\,D_{t}\frac{\mathfrak{v}_{1}}{2}=D\sigma_{1}$

then it has the form

(57) $\displaystyle u_{t}$ $\displaystyle=v,$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+g_{u_{2}}u_{3}^{2}+2g_{v}v_{1}u_{3}+2\left(g_{u_{1}}u_{2}+g_{u}u_{1}\right)u_{3}+f_{2}v_{1}^{2}+f_{1}v_{1}+f_{0}.$

where $f_{i}=f_{i}(u,u_{1},u_{2},v)$.

###### Proof.

Condition (55) amounts to

$f_{u_{3}}(u,u_{1},u_{2},u_{3},v,v_{1})=2Dg(u,u_{1},u_{2},v)=2g_{u_{2}}u_{3}+2g_{v}v_{1}+2g_{u_{1}}u_{2}+2g_{u}u_{1}$

and, therefore

$f=g_{u_{2}}u_{3}^{2}+2g_{v}v_{1}u_{3}+2\left(g_{u_{1}}u_{2}+g_{u}u_{1}\right)u_{3}+\tilde{f}(u,u_{1},u_{2},v,v_{1}).$

Now condition (56) becomes equivalent to

(58)
$D_{t}(2g_{v}u_{3}+\tilde{f}_{v_{1}})=D\sigma_{1}(u,u_{1},u_{2},u_{3},u_{4},v,v_{1},v_{2}),$

where $\sigma_{1}$ is an unknown function. By subtracting an appropriate total
$x$-derivative from both sides, we can lower the differential order of the
left-hand side until the dependence of the right-hand side on higher
derivatives becomes nonlinear. In particular,

$D_{t}(2g_{v}u_{3}+\tilde{f}_{v_{1}})-D(u_{4}\tilde{f}_{v_{1},v_{1}}+2g_{v}v_{2})=-\tilde{f}_{v_{1}v_{1}v_{1}}u_{4}v_{2}+Au_{4}+Bv_{2}+C$

where $A$, $B$, $C$ are functions of $(u,u_{1},u_{2},u_{3},v,v_{1})$. Since
the rhs has to be a total $x$-derivative, it should be
$\tilde{f}_{v_{1}v_{1}v_{1}}=0$. This proves the lemma and illustrates the
general procedure to solve the system of integrability conditions. ∎

######  Remark 5.

The standard way to eliminate $\sigma_{1}$ from (58) is applying Euler
operators (cf. e.g. [3]) to both sides of the equation. We would obtain

$0=\frac{\delta}{\delta
u}\Big{(}D_{t}(2g_{v}u_{3}+\tilde{f}_{v_{1}})\Big{)}=-\frac{1}{2}\tilde{f}_{v_{1}v_{1}v_{1}}v_{6}-\frac{5}{2}\tilde{f}_{v_{1}v_{1}v_{1}v_{1}}v_{2}v_{5}\\\
+\left(g_{u_{2}}\tilde{f}_{v_{1}v_{1}v_{1}}-\frac{5}{2}\tilde{f}_{u_{2}v_{1}v_{1}v_{1}}\right)u_{3}v_{5}+\left(g_{u_{2}}\tilde{f}_{v_{1}v_{1}v_{1}}-\frac{3}{2}\tilde{f}_{u_{2}v_{1}v_{1}v_{1}}\right)v_{2}u_{6}\\\
-\left(6g_{vv}g_{u_{2}}-3\tilde{f}_{u_{2}v_{1}v_{1}}g_{u_{2}}-3g_{u_{2}vv}+\frac{3}{2}\tilde{f}_{u_{2}u_{2}v_{1}v_{1}}\right)u_{3}u_{6}+\cdots\\\
\Rightarrow\tilde{f}_{v_{1}v_{1}v_{1}}=0.$

We prefer the procedure based on the refinement of $\sigma_{1}$ for two
reasons. First, to use the condition $C_{3}$ we need to know at least the
dependence of $\sigma_{1}$ on higher derivatives. Elimination of $\sigma_{1}$
won’t help here. Second, the direct computation of the variational derivative
of long expressions is very demanding on computational resources.

Continuing the classification of integrable systems (72) requires considering
a function $g$ of different types. In this paper we investigate the cases
$g=g(u)$ and $g=g(u_{1})$ (it can be proved that the case $g=g(u,u_{1})$ can
be reduced to these two subcases). It turns out that, even in these simple
cases, we arrive at non-trivial lists of integrable equations.

######  Remark 6.

The class of systems (72) is invariant with respect to point transformations

(59) $u\to\varphi(p),\qquad v\to\varphi^{\prime}(p)q.$

A classification of integrable systems (72) should be performed up to
transformations (59). It is easy to verify that in the case $g=g(u,u_{1})$ the
function $g$ transforms as:

$g(u,u_{1})\to
2\log\varphi^{\prime}(p)+g(\varphi(p),\varphi^{\prime}(p)p_{1}).$

##  6\. Classification of the case $g=g(u)$

Using a transformation (59), we reduce $g$ to zero. Then the system (72)
becomes

(60) $\displaystyle u_{t}$ $\displaystyle=v,$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+f_{2}v_{1}^{2}+f_{1}v_{1}+f_{0},$

with $f_{i}=f_{i}(u,u_{1},u_{2},v)$, $i=0,1,2$.

######  Theorem 2.

If a nonlinear equation of the form (60) satisfies the integrability
conditions $C_{0}$–$C_{8}$ then it can be reduced, using transformations
$t\to\alpha t$, $x\to\beta x+\gamma t$, $u\to\delta u+\kappa x+\lambda t+\xi
x^{2}+\chi t^{2}$ to one of the following:

(61) $\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+u_{2}^{2}\end{aligned}\right.$ (62)
$\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+2u_{1}u_{2}\end{aligned}\right.$ (63)
$\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+2u_{2}v+2u_{1}^{2}u_{2}\end{aligned}\right.$ (64)
$\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+4u_{1}v_{1}+2u_{2}v-6u_{1}^{2}u_{2}.\end{aligned}\right.$

In order to prove Theorem 2, we first perform a preliminary classification
producing an ansatz for integrable systems (60) where the only undetermined
expressions are constants. In a second step we continue applying integrability
conditions until obtaining a final list of equations (61)-(64) that pass a
large (but finite) number of conditions. The third step is to directly prove
the integrability of the listed equations by finding a recursion operator for
each of them (see Section 8).

######  Lemma 2.

All integrable systems of the form (60) can be written as

(65) $\displaystyle u_{t}$ $\displaystyle=v$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+(\alpha_{0}+\alpha_{1}u_{1})v_{1}+2\left(\beta-\frac{1}{4}\alpha_{1}\right)u_{2}v+(\delta_{0}+\delta_{1}u)u_{2}^{2}$
$\displaystyle\qquad{}+\frac{1}{8}\left(16\delta_{1}-\alpha_{1}^{2}-12\alpha_{1}\beta+16\beta^{2}\right)u_{1}^{2}u_{2}+(\kappa_{0}+\kappa_{1}u_{1})u_{2}+\tau_{0}$

where $\alpha_{0}$, $\alpha_{1}$, $\beta$, $\delta_{0}$, $\delta_{1}$,
$\kappa_{0}$, $\kappa_{1}$, $\tau_{0}$ are constants.

###### Scheme of the proof of Lemma 2.

It follows from $C_{1}$, $C_{2}$ and $C_{3}$ that

(66) $\displaystyle f_{2}(u,u_{1},u_{2},v)=0,$ (67) $\displaystyle
f_{1}(u,u_{1},u_{2},v)=a(u,u_{1})+f_{10}(u,u_{1})v+f_{11}(u,u_{1})u_{2}.$

Now condition $C_{2}$ allows to write $f_{0}$ in terms of other functions as

(68)
$f_{0}(u,u_{1},u_{2},v)=\tilde{f}_{0}\left(u,u_{1},u_{2}\right)+\tfrac{1}{2}\left[(f_{10})_{u_{1}}u_{2}+(f_{10})_{u}u_{1}\right]v^{2}+2(h_{u_{1}}u_{2}+h_{u}u_{1})v$

where $h$ is a function of $u$ and $u_{1}$. Condition $C_{3}$ implies that the
function $f_{0}(u,u_{1},u_{2})$ can be written as

$\tilde{f}_{0}(u,u_{1},u_{2})=f_{00}+u_{2}f_{01}+u_{2}^{2}f_{02}-\tfrac{1}{8}u_{2}^{3}f_{10}^{2}-\tfrac{1}{24}u_{2}^{3}f_{11}^{2}$

where $f_{00}$, $f_{01}$ and $f_{02}$ are functions of $u$ and $u_{1}$.
Condition $C_{4}$ sets $f_{11}(u,u_{1})=\bar{f}_{11}(u)$ and from $C_{3}$
arises the following condition:

$f_{10}(u,u_{1})\bar{f}_{11}(u)=0.$

Setting $\bar{f}_{11}(u)=0$ implies, using $C_{6}$, that $f_{10}(u,u_{1})=0$.
So we have that anyway $f_{10}(u,u_{1})=0$. Then $C_{1}$ is equivalent to

$a(u,u_{1})=\alpha_{0}+a_{1}(u)u_{1}+\bar{f}^{\prime}_{11}(u)u_{1}^{2}.$

Using in a straightforward way conditions $C_{3}$, $C_{4}$, $C_{5}$ and
$C_{6}$ one can find that

$\displaystyle\bar{f}_{11}(u)=0,\qquad
a(u,u_{1})=\alpha_{0}+\alpha_{1}\,u_{1},$ $\displaystyle
h(u,u_{1})=\left(\beta-\tfrac{1}{4}\alpha_{1}\right)\,u_{1},$ $\displaystyle
f_{00}(u,u_{1})=\bar{f}_{00}(u)$
$\displaystyle\begin{split}f_{01}(u,u_{1})&=-\frac{1}{8}u_{1}^{2}\left(-16\delta_{1}+\alpha_{1}^{2}+12\alpha_{1}\,\beta-16\beta^{2}\right)+\kappa_{0}+\kappa_{1}u_{1},\end{split}$
$\displaystyle f_{02}(u,u_{1})=\delta_{0}+\delta_{1}u,$

where $\beta$, $\alpha_{0}$, $\alpha_{1}$, $\delta_{0}$, $\delta_{1}$,
$\kappa_{0}$, $\kappa_{1}$ are arbitrary constants, and the branching
condition

$\beta\cdot\bar{f}_{00}(u)=0.$

At this moment the system has the form

$\displaystyle u_{t}$ $\displaystyle=v$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+(\alpha_{0}+\alpha_{1}u_{1})v_{1}+(\delta_{0}+\delta_{1}u)u_{2}^{2}+2\left(\beta-\frac{1}{4}\alpha_{1}\right)u_{2}v$
$\displaystyle\qquad{}+\frac{1}{8}\left(16\delta_{1}-\alpha_{1}^{2}-12\alpha_{1}\beta+16\beta^{2}\right)u_{1}^{2}u_{2}+(\kappa_{0}+\kappa_{1}u_{1})u_{2}+\bar{f}_{00}(u).$

If we suppose that $\bar{f}_{00}(u)\neq 0$, the conditions up to $C_{9}$
require that $\bar{f}_{00}(u)$ is constant, so the lemma is proved. ∎

###### Scheme of the proof of Theorem 2.

Imposing the integrability conditions up to $C_{9}$ to a system of type (65)
yields a system of algebraic equations over the constants $\beta$,
$\alpha_{0}$, $\alpha_{1}$, $\delta_{0}$, $\delta_{1}$, $\kappa_{0}$,
$\kappa_{1}$ equivalent to

$\displaystyle\delta_{1}=0,$
$\displaystyle\alpha_{1}\delta_{0}=\kappa_{1}\delta_{0}=\beta_{1}\delta_{0}=0,$
$\displaystyle\alpha_{1}\tau_{0}=\kappa_{1}\tau_{0}=\beta_{1}\tau_{0}=0,$
$\displaystyle\alpha_{1}\left(\alpha_{1}-2\beta_{1}\right)=0,$
$\displaystyle\alpha_{1}\left(3\alpha_{0}\alpha_{1}+4\kappa_{1}\right)=0.$

Solving the previous system leads to the following three subclasses of
systems:

$\displaystyle\hskip 15.0pt\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\alpha_{0}v_{1}+\delta_{1}u_{2}^{2}+\kappa_{0}u_{2}+\tau_{0},\end{aligned}\right.$
$\displaystyle\hskip 15.0pt\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\alpha_{0}v_{1}+2\beta
vu_{2}+2\beta^{2}u_{1}^{2}u_{2}+\kappa_{1}u_{1}u_{2}+\kappa_{0}u_{2},\end{aligned}\right.$
$\displaystyle\hskip 15.0pt\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\alpha_{1}u_{1}v_{1}+\alpha_{0}v_{1}+\frac{1}{2}\alpha_{1}vu_{2}-\frac{3}{8}\alpha_{1}^{2}u_{1}^{2}u_{2}-\frac{3}{4}\alpha_{0}\alpha_{1}u_{1}u_{2}+\kappa_{0}u_{2}.\end{aligned}\right.$

A Galilean transformation $x\to x+\alpha t$ can set $\alpha_{0}=0$, and a
further change $u\to u+\xi x^{2}$ allows to eliminate $\kappa_{0}$. The
constant $\tau_{0}$ can be cancelled by the change $u\to u+\tau_{0}\,t^{2}/2$,
$v\to v+\tau_{0}\,t$ so the systems become

(69) $\displaystyle\hskip 15.0pt\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\delta_{1}u_{2}^{2},\end{aligned}\right.$ (70)
$\displaystyle\hskip 15.0pt\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+2\beta
vu_{2}+2\beta^{2}u_{1}^{2}u_{2}+\kappa_{1}u_{1}u_{2},\end{aligned}\right.$
(71) $\displaystyle\hskip 15.0pt\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\alpha_{1}u_{1}v_{1}+\frac{1}{2}\alpha_{1}vu_{2}-\frac{3}{8}\alpha_{1}^{2}u_{1}^{2}u_{2}.\end{aligned}\right.$

The system (69) can be scaled to (61). The system (70) contains two normalized
cases: when $\beta=0$, $u\to\delta u+\kappa x$ leads to (62), while when
$\beta\neq 0$ a transformation $x\to\gamma x$, $t\to\gamma^{2}t$, $u\to\delta
u+\kappa x+\lambda t$ leads to (63). The system (71) can be reduced to (64)
with a transformation $x\to x+\tau t$, $u\to\delta u+\lambda t$. ∎

##  7\. Classification of the case $g=g(u_{1})$

When $g=g(u_{1})$ system (72) becomes

(72) $\displaystyle u_{t}$ $\displaystyle=v,$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+2g^{\prime}(u_{1})u_{2}u_{3}+f_{2}v_{1}^{2}+f_{1}v_{1}+f_{0}.$

where $f_{i}=f_{i}(u,u_{1},u_{2},v)$. In this section we consider that
$g^{\prime}(u_{1})\neq 0$, to avoid falling into the previous case.

######  Theorem 3.

If a nonlinear equation of the form (72) satisfies the integrability
conditions $C_{0}$–$C_{10}$ then it can be reduced, using transformations
$t\to\alpha t$, $x\to\beta x+\gamma t$, $u\to\delta u+\kappa x+\lambda t$ to
one of the following

(73) $\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}-\frac{4}{u_{1}}u_{2}u_{3}+\frac{1}{u_{1}^{2}}u_{2}v^{2}+\frac{3}{u_{1}^{2}}u_{2}^{3}\end{aligned}\right.$
(74) $\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}-\frac{4}{u_{1}}u_{2}u_{3}+\frac{1}{u_{1}^{2}}u_{2}v^{2}+\frac{3}{u_{1}^{2}}u_{2}^{3}+\left(\epsilon
u_{1}^{2}+\frac{1}{u_{1}}\right)u_{2}\end{aligned}\right.$ (75)
$\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}-\frac{4}{u_{1}}u_{2}u_{3}+\frac{1}{u_{1}^{2}}u_{2}v^{2}-\frac{2}{u_{1}^{2}}u_{2}v+\frac{3}{u_{1}^{2}}u_{2}^{3}+\left(\frac{\epsilon}{u_{1}}+\frac{1}{u_{1}^{2}}\right)u_{2}\end{aligned}\right.$
(76) $\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}-\frac{4}{u_{1}}u_{2}u_{3}+\frac{1}{u_{1}^{2}}u_{2}v^{2}-\frac{2}{u_{1}^{2}}u_{2}v+\frac{3}{u_{1}^{2}}u_{2}^{3}+\left(u_{1}^{2}+\frac{c}{u_{1}}+\frac{1}{u_{1}^{2}}\right)u_{2}\end{aligned}\right.$
(77) $\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\frac{4}{u_{1}}vv_{1}-\frac{4}{u_{1}}u_{2}u_{3}-\frac{3}{u_{1}^{2}}u_{2}v^{2}+\frac{3}{u_{1}^{2}}u_{2}^{3}\end{aligned}\right.$
(78) $\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\frac{4}{u_{1}}vv_{1}+\frac{4}{u_{1}}v_{1}-\frac{4}{u_{1}}u_{2}u_{3}-\frac{3}{u_{1}^{2}}u_{2}v^{2}-\frac{2}{u_{1}^{2}}u_{2}v+\frac{3}{u_{1}^{2}}u_{2}^{3}\\\
&\quad{}+\left(\epsilon
u_{1}^{2}-\frac{1}{3u_{1}^{2}}\right)u_{2}\end{aligned}\right.$ (79)
$\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\frac{4u_{1}}{u_{1}^{2}+1}vv_{1}-\frac{4u_{1}}{u_{1}^{2}+1}u_{2}u_{3}-\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}v^{2}+\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}^{3}\\\
&\quad{}+\left(u_{1}^{2}+1\right)s(u)u_{2}+\frac{1}{6}s^{\prime}(u)\left(u_{1}^{2}+1\right)\left(3u_{1}^{2}-1\right)\\\
&\text{with
}s^{\prime\prime\prime}(u)+4s(u)s^{\prime}(u)=0\end{aligned}\right.$ (80)
$\displaystyle\left\\{\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\frac{4u_{1}}{u_{1}^{2}+1}vv_{1}+\frac{2u_{1}}{u_{1}^{2}+1}v_{1}-\frac{4u_{1}}{u_{1}^{2}+1}u_{2}u_{3}-\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}v^{2}\\\
&\quad{}-\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}v+\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}^{3}-\frac{3u_{1}^{2}-1}{4\left(u_{1}^{2}+1\right)^{2}}u_{2}+c\left(u_{1}^{2}+1\right)u_{2}\end{aligned}\right.$

where $\epsilon$ is equal to 0 or 1 and $c$ is an arbitrary constant.

###### Scheme of the proof of Theorem 3.

Formulas (66) and (67) remain valid in this case. Equation (68) becomes

$f_{0}(u,u_{1},u_{2},v)=\tilde{f}_{0}(u,u_{1},u_{2})+\left[\left(\tfrac{1}{2}g_{u_{1}u_{1}}+\tfrac{1}{4}g_{u_{1}}f_{10}+\tfrac{1}{2}(f_{10})_{u_{1}}\right)u_{2}+\tfrac{1}{2}(f_{10})_{u}u_{1}\right]v^{2}\\\
{}+\left[\tfrac{1}{2}g_{u_{1}}f_{11}u_{2}^{2}+\tfrac{1}{2}g_{u_{1}}au_{2}+2h_{u_{1}}u_{2}+2h_{u}u_{1}\right]v$

and now

$\tilde{f}_{0}(u,u_{1},u_{2})=f_{00}+u_{2}f_{01}+u_{2}^{2}f_{02}-\tfrac{1}{8}f_{10}^{2}u_{2}^{3}-\tfrac{1}{24}f_{11}^{2}u_{2}^{3}\\\
{}-\tfrac{1}{4}g_{u_{1}}f_{10}u_{2}^{3}+\tfrac{1}{2}\left(g_{u_{1}}\right)^{2}u_{2}^{3}+\tfrac{1}{2}g_{u_{1}u_{1}}u_{2}^{3}.$

When $g_{u_{1}}(u_{1})\neq 0$, as is the case, $C_{4}$ implies that
$f_{11}(u,u_{1})=0$. Searching up to $C_{8}$ we find the following two branch
conditions:

$\displaystyle\left[f_{10}+2g_{u_{1}}\right]\left[(f_{10})_{u_{1}}-g_{u_{1}}f_{10}\right]=0,$
$\displaystyle\left[f_{10}+2g_{u_{1}}\right]\left[2(f_{10})_{u_{1}}+f_{10}^{2}\right]=0.$

###  7.1. Case $f_{10}(u,u_{1})+2g^{\prime}(u_{1})=0$

Condition $C_{1}$ implies that $f_{02}(u,u_{1})=0$ and that the function
$g(u_{1})$ must satisfy the equation

$g^{\prime\prime\prime}-3g^{\prime}g^{\prime\prime}+(g^{\prime})^{2}=0$

which means that, without loss of generality, we can set

$g(u_{1})=\log{\frac{1}{\gamma_{0}+(u_{1}+\gamma_{1})^{2}}}.$

Conditions $C_{3}$, $C_{4}$ imply that we can set

$h(u,u_{1})=\frac{1}{4}a(u,u_{1}).$

Now two subcases must be considered, $\gamma_{0}=0$ and $\gamma_{0}\neq 0$.

####  7.1.1. Subcase $\gamma_{0}=0$

From $C_{1}$

$a(u,u_{1})=\alpha_{0}+\frac{b(u)}{u_{1}+\gamma_{1}}$

and from $C_{5}$

$f_{01}(u,u_{1})=c(u)(u_{1}+\gamma_{1})^{2}-\frac{\alpha_{0}}{4}\frac{b(u)}{u_{1}+\gamma_{1}}-\frac{3}{16}\frac{b(u)^{2}}{(u_{1}+\gamma_{1})^{2}}.$

Now $C_{1}$ imposes that

$f_{00}\left(u,u_{1}\right)=-\frac{1}{6}c^{\prime}(u)u_{1}\left(\gamma_{1}+u_{1}\right)\left(5\gamma_{1}^{2}-6\gamma_{1}u_{1}-3u_{1}^{2}\right)\\\
{}+\frac{1}{4}\alpha_{0}b^{\prime}(u)u_{1}+\frac{3b(u)b^{\prime}(u)u_{1}}{16\left(\gamma_{1}+u_{1}\right)}.$

Again a branch appears, with $\gamma_{1}=0$ or $\gamma_{1}\neq 0$.

When $\gamma_{1}=0$ all conditions up to $C_{10}$ are solved only if

(81) $b(u)c^{\prime}(u)+2b^{\prime}(u)c(u)+2b^{\prime\prime\prime}(u)=0.$

The system has the form

$\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\left(\alpha_{0}+\frac{b(u)}{u_{1}}+\frac{4v}{u_{1}}\right)v_{1}-\frac{4u_{2}}{u_{1}}u_{3}-\frac{3u_{2}}{u_{1}^{2}}v^{2}\\\
&\quad{}-\left(\frac{\alpha_{0}}{u_{1}}+\frac{3}{2}\frac{b(u)}{u_{1}^{2}}\right)u_{2}v+\frac{1}{2}b^{\prime}(u)v+\frac{3u_{2}^{3}}{u_{1}^{2}}+\left(c(u)u_{1}^{2}-\frac{\alpha_{0}b(u)}{4u_{1}}-\frac{3b(u)^{2}}{16u_{1}^{2}}\right)u_{2}\\\
&\quad{}+\frac{1}{16}\left(8c(u)_{u}u_{1}^{4}+4\alpha_{0}b^{\prime}(u)u_{1}+3b(u)b^{\prime}(u)\right)\end{aligned}$

A Galilean transformation eliminates $\alpha_{0}$. If $b(u)=0$, a
transformation $\varphi(u)\to u$ with

$2\frac{\varphi^{\prime\prime\prime}(u)}{\varphi^{\prime}(u)}-3\frac{\varphi^{\prime\prime}(u)^{2}}{\varphi^{\prime}(u)^{2}}=c(u),$

produces equation (77). If $b(u)\neq 0$ a transformation $\varphi(u)\to u$
with $\varphi^{\prime}(u)=1/b(u)$, condition (81) and appropriate scalings
lead to equation (78).

When $\gamma_{1}\neq 0$ we have that $C_{1}$ implies that $b(u)$ and $c(u)$
are constant, which means that the expression of the system has no function
explicitly dependent on $u$. Then we can perform a transformation $u\to
u-\gamma_{1}x$ and set $\gamma_{1}=0$, so we fall into the case studied
immediately above.

####  7.1.2. Subcase $\gamma_{0}\neq 0$

We normalize $\gamma_{0}=1$. From $C_{1}$ and $C_{3}$

$a(u,u_{1})=\alpha_{0}\frac{1-(u_{1}+\gamma_{1})^{2}}{1+(u_{1}+\gamma_{1})^{2}}+\alpha_{1}\frac{2\left(\gamma_{1}+u_{1}\right)}{1+(u_{1}+\gamma_{1})^{2}}.$

Now from $C_{1}$ we see that

$\displaystyle
f_{01}(u,u_{1})=\frac{1+(u_{1}+\gamma_{1})^{2}}{2(u_{1}+\gamma_{1})}\left(\sigma_{1}(u,u_{1})\right)_{u_{1}},$
$\displaystyle
f_{00}(u,u_{1})=\frac{1+(u_{1}+\gamma_{1})^{2}}{2(u_{1}+\gamma_{1})}\left(\sigma_{1}(u,u_{1})\right)_{u}u_{1}$

being $\sigma_{1}(u,u_{1})$ an arbitrary function. From $C_{3}$ a branching
condition arises

$\alpha_{1}\left(\sigma_{1}(u,u_{1})\right)_{u}=0.$

When $\alpha_{1}=0$, from conditions $C_{3}$ and $C_{5}$ it follows that

$\sigma_{1}(u,u_{1})=-\frac{\alpha_{0}^{2}\left(1+3\left(u_{1}+\gamma_{1}\right)^{2}\right)}{4\left(1+\left(u_{1}+\gamma_{1}\right)^{2}\right)^{2}}+\frac{1}{2}\sigma(u)\left(\left(u_{1}+\gamma_{1}\right)^{2}-\frac{1}{3}\right)+\tau_{0}$

and another branching condition

$\gamma_{1}\cdot\sigma^{\prime}(u)=0.$

If $\gamma_{1}=0$, $C_{5}$ implies that

(82) $\sigma^{\prime\prime\prime}+2\sigma\sigma^{\prime}=0$

and all the conditions up to $C_{8}$ are satisfied. The resulting equation is

$\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\frac{4u_{1}}{u_{1}^{2}+1}vv_{1}-\frac{\alpha_{0}\left(u_{1}^{2}-1\right)}{u_{1}^{2}+1}v_{1}-\frac{4u_{1}}{u_{1}^{2}+1}u_{2}u_{3}\\\
&\quad{}-\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}v^{2}+\frac{\alpha_{0}u_{1}\left(u_{1}^{2}-3\right)}{\left(u_{1}^{2}+1\right)^{2}}u_{2}v+\frac{\left(3u_{1}^{2}-1\right)u_{2}^{3}}{\left(u_{1}^{2}+1\right)^{2}}\\\
&\quad{}+\frac{\alpha_{0}^{2}(3u_{1}^{2}-1)}{4\left(u_{1}^{2}+1\right)^{2}}u_{2}+\frac{1}{2}\left(u_{1}^{2}+1\right)\sigma(u)u_{2}\\\
&\quad{}+\frac{1}{12}\sigma^{\prime}(u)\left(u_{1}^{2}+1\right)\left(3u_{1}^{2}-1\right)\end{aligned}$

that after Galilei and $\sigma(u)\to 2s(u)$ becomes

$\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\frac{4u_{1}}{u_{1}^{2}+1}vv_{1}-\frac{4u_{1}}{u_{1}^{2}+1}u_{2}u_{3}-\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}v^{2}+\frac{\left(3u_{1}^{2}-1\right)u_{2}^{3}}{\left(u_{1}^{2}+1\right)^{2}}\\\
&\quad{}+\left(u_{1}^{2}+1\right)s(u)u_{2}+\frac{1}{6}s^{\prime}(u)\left(u_{1}^{2}+1\right)\left(3u_{1}^{2}-1\right)\end{aligned}$

with $\sigma^{\prime\prime\prime}+4\sigma\sigma^{\prime}=0$. This is eq. (79).

If $\gamma_{1}\neq 0$, we have that $\sigma(u)=\text{constant}$ and no more
restrictions up to level $C_{8}$. The resulting equation, with Galilei,
scaling and $u\to u+\alpha x$ transformations, can be reduced to (79) with
constant $s(u)$.

When $\alpha_{1}\neq 0$, we have that $\sigma_{1}(u,u_{1})=\sigma(u_{1})$ and
Galilei, scaling and $u\to u+\alpha x$ transformations lead to

$\begin{aligned} u_{t}&=v\\\
v_{t}&=u_{4}+\frac{4u_{1}}{u_{1}^{2}+1}vv_{1}+\frac{2\alpha_{1}u_{1}}{u_{1}^{2}+1}v_{1}-\frac{4u_{1}u_{2}}{u_{1}^{2}+1}u_{3}-\frac{\left(3u_{1}^{2}-1\right)u_{2}}{\left(u_{1}^{2}+1\right)^{2}}v^{2}\\\
&\quad{}-\frac{\alpha_{1}\left(3u_{1}^{2}-1\right)u_{2}}{\left(u_{1}^{2}+1\right)^{2}}v+\frac{3u_{1}^{2}-1}{\left(u_{1}^{2}+1\right)^{2}}u_{2}^{3}-\frac{\alpha_{1}^{2}\left(3u_{1}^{2}-1\right)}{4\left(u_{1}^{2}+1\right)^{2}}u_{2}+\sigma_{0}\left(u_{1}^{2}+1\right)u_{2}\end{aligned}$

Only scalings are available to simplify the equation now, and we obtain a
subcase of (79) and eq. (80).

###  7.2. Case $f_{10}(u,u_{1})+2g_{u_{1}}(u_{1})\neq 0$

We have that $f_{10}(u,u_{1})=0$. Conditions up to $C_{6}$ imply that
$f_{02}(u,u_{1})=0$ and that $g(u_{1})$ must satisfy

$2g^{\prime\prime}-(g^{\prime})^{2}=0$

i.e. we can take

$g(u_{1})=-2\log(u_{1}+\gamma_{1}).$

Conditions up to $C_{7}$ imply that

$\displaystyle a(u,u_{1})=\alpha_{0}\text{ constant},$ $\displaystyle
h(u,u_{1})=\frac{d(u)}{u_{1}+\gamma_{1}}.$

The form of the system at this stage is

$\displaystyle u_{t}$ $\displaystyle=v$ $\displaystyle v_{t}$
$\displaystyle=u_{4}+\alpha_{0}v_{1}-\frac{4u_{3}u_{2}}{u_{1}+\gamma_{1}}+\frac{v^{2}u_{2}}{\left(u_{1}+\gamma_{1}\right)^{2}}-\frac{\alpha_{0}vu_{2}}{u_{1}+\gamma_{1}}-\frac{2d(u)vu_{2}}{\left(u_{1}+\gamma_{1}\right)^{2}}$
$\displaystyle\quad{}+\frac{3u_{2}^{3}}{\left(u_{1}+\gamma_{1}\right)^{2}}+\frac{2d^{\prime}(u)vu_{1}}{u_{1}+\gamma_{1}}+u_{2}f_{01}\left(u,u_{1}\right)+f_{00}\left(u,u_{1}\right).$

Condition $C_{3}$ implies that

$f_{01}\left(u,u_{1}\right)_{u_{1}u_{1}}-\frac{2f_{01}\left(u,u_{1}\right)}{\left(u_{1}+\gamma_{1}\right)^{2}}-\frac{4d(u)^{2}}{\left(u_{1}+\gamma_{1}\right)^{4}}=0$

so we write

$f_{01}(u,u_{1})=c(u)\left(u_{1}+\gamma_{1}\right)^{2}+\frac{b(u)}{u_{1}+\gamma_{1}}+\frac{d(u)^{2}}{\left(u_{1}+\gamma_{1}\right)^{2}}.$

Besides we obtain

$f_{00}(u,u_{1})=\frac{1}{2}c^{\prime}(u)\left(u_{1}-\gamma_{1}\right)\left(u_{1}+\gamma_{1}\right)^{3}-\frac{1}{2}b^{\prime}(u)\left(2u_{1}+\gamma_{1}\right)-d(u)d^{\prime}(u)$

and the condition

$\gamma_{1}\cdot d^{\prime}(u)=0.$

Supposing that $\gamma_{1}\neq 0$ implies that $d(u)$ is constant and, by
means of $C_{4}$, that also $b(u)$ and $c(u)$ are constant. That makes the rhs
of the system not explicitly dependent on $u$, and thus a transformation $u\to
u-\gamma_{1}x$ can set $\gamma_{1}$ to zero, without introducing explicit
dependencies on $x$. We can thus set $\gamma_{1}=0$, yielding

$\displaystyle
f_{01}(u,u_{1})=c(u)u_{1}^{2}+\frac{b(u)}{u_{1}}+\frac{d(u)^{2}}{u_{1}^{2}},$
$\displaystyle
f_{00}(u,u_{1})=\frac{1}{2}c^{\prime}(u)u_{1}^{4}-b^{\prime}(u)u_{1}-d(u)d^{\prime}(u).$

Condition $C_{3}$ implies a relation

$d(u)b^{\prime}(u)-d^{\prime}(u)b(u)=0$

opening two subcases, $d(u)=0$ and $d(u)\neq 0$.

####  7.2.1. Case $d(u)\neq 0$

We have $b(u)=\beta d(u)$ and

$d(u)c^{\prime}(u)+2d^{\prime}(u)c(u)+2d^{\prime\prime\prime}(u)=0$

so we can write

$c(u)=\frac{\gamma}{d^{2}(u)}-\frac{2d(u)d^{\prime\prime}(u)-(d^{\prime}(u))^{2}}{d^{2}(u)}$

We obtain thus a system whose rhs depends on one arbitrary function $d(u)$ and
some arbitrary constants, and that passes all integrability conditions up to
$C_{10}$. A transformation $\varphi(u)\to u$ with $\varphi^{\prime}(u)=1/d(u)$
can set $d(u)=1$ yielding, after some scalings and Galilei, the system

(83) $\displaystyle u_{t}$ $\displaystyle=v$ $\displaystyle v_{t}$
$\displaystyle=u_{4}-\frac{4u_{3}u_{2}}{u_{1}}+\frac{v^{2}u_{2}}{u_{1}^{2}}-\frac{2vu_{2}}{u_{1}^{2}}+\frac{3u_{2}^{3}}{u_{1}^{2}}+\left(\kappa_{0}u_{1}^{2}+\frac{\kappa_{1}}{u_{1}}+\frac{1}{u_{1}^{2}}\right)u_{2}$

Different normalizations and scalings produce equations (75) and (76).

####  7.2.2. Case $d(u)=0$

Condition $C_{5}$ yields again the equation

$b(u)c^{\prime}(u)+2b^{\prime}(u)c(u)+2b^{\prime\prime\prime}(u)=0$

which can be solved, if $b(u)\neq 0$, in $c(u)$. The resulting system can be
again transformed through $\varphi(u)\to u$ with $\varphi^{\prime}(u)=1/b(u)$
into (74). If $b(u)=0$, a Galilean transformation and $\varphi(u)\to u$ with
$2\frac{\varphi^{\prime\prime\prime}(u)}{\varphi^{\prime}(u)}-3\frac{\varphi^{\prime\prime}(u)^{2}}{\varphi^{\prime}(u)^{2}}=c(u)$
allows to get the system (73). ∎

######  Remark 7.

In [18] a list of integrable Lagrangian systems with Lagrangian density
$\mathscr{L}=\frac{1}{2}L_{2}(u,u_{1},u_{2})u_{t}^{2}+L_{1}(u,u_{1},u_{2})u_{t}+L_{0}(u,u_{1},u_{2})$
was found. Comparing with that list, we find that Eq. (63) is equivalent to
the Lagrangian system L3, Eq. (64) to L2 and Eq. (79) to L8.

##  8\. Recursion, symplectic and Hamiltonian operators for the equations in
Theorem 2

The systems found in Theorems 2 and 3 satisfy only a finite number of
necessary conditions for integrability (46). But integrability requires the
satisfaction of an infinite number of integrability conditions. To prove
integrability we will determine explicit recursion and symplectic operators
for each of the listed systems.

Recall that there are the following algebraic relations between operators
[19]. For any recursion $\mathbf{R}$, symplectic $\mathbf{S}$ and Hamiltonian
$\mathbf{H}$ operators, the product $\mathbf{S}\mathbf{R}$ is a symplectic
operator, the product $\mathbf{R}\mathbf{H}$ is a Hamiltonian operator and the
inverse $\mathbf{S}^{-1}$ is a Hamiltonian operator (and viceversa, the
inverse $\mathbf{H}^{-1}$ is a symplectic operator).

Weakly nonlocal operators. Almost all known recursion, symplectic and
Hamiltonian operators for various integrable systems can be written in quasi-
local (or weakly nonlocal) form.

######  Definition 4.

A matrix pseudo-differential operator $\mathbf{A}$ is said to be written in a
quasi-local form if

(84)
$\mathbf{A}=\mathbf{D}+\sum_{i=1}^{k}\mathbf{B}_{i}D^{-1}\cdot\mathbf{C}_{i}^{+}$

where $\mathbf{D}$ is a differential operator and $\mathbf{B}_{i}$,
$\mathbf{C}_{i}$ are rectangular matrices of differential functions of
compatible size.

The matrices $\mathbf{B}_{i}$, $\mathbf{C}_{i}$ are related [14, 15, 17, 16]
to local symmetries and cosymmetries of the corresponding integrable system.
Recall [3] that symmetries and cosymmetries are vector-columns $\bf b$ and
$\bf c$, depending on jet variables, satisfying respectively the equations

$\mathbf{b}_{t}=\mathbf{\Phi}_{*}(\mathbf{b}),\qquad\mathbf{c}_{t}=-\mathbf{\Phi}_{*}^{+}(\mathbf{c}).$

Here the $t$-derivative is calculated in virtue of the system under
consideration. A special type of cosymmetry is the variational derivative of a
conserved density $\rho$, that is of the form

$\mathbf{c}=\frac{\delta\rho}{\delta{\mathbf{u}}}=\begin{bmatrix}\delta\rho\big{/}\delta
u\\\ \delta\rho/\delta v\end{bmatrix},$

where

$\frac{\delta\rho}{\delta u}=\sum_{i\geq
0}D^{i}\left(\frac{\partial\rho}{\partial
u_{i}}\right),\quad\frac{\delta\rho}{\delta v}=\sum_{i\geq
0}D^{i}\left(\frac{\partial\rho}{\partial v_{i}}\right).$

In the case of recursion operators, in formula (84) the columns of matrices
$\mathbf{B}_{i}$ and $\mathbf{C}_{i}$ are correspondingly symmetries and
cosymmetries. For symplectic operators both $\mathbf{B}_{i}$ and
$\mathbf{C}_{i}$ consist of cosymmetries, and for Hamiltonian operators both
$\mathbf{B}_{i}$ and $\mathbf{C}_{i}$ are formed by symmetries.

Quasi-local anzats are practical because finding symmetries and cosymmetries
of a system is relatively simple. If the rhs of a system as well as its
symmetries and cosymmetries are homogeneous polynomials, computations become
almost trivial.

In the simplest version of the ansatz (84) the matrices $\mathbf{B}_{i}$ and
$\mathbf{C}_{i}$ are column-vectors. We will call such operators vector quasi-
local operator. In particular, vector quasi-local recursion operator
$\mathbf{R}$ has the form

(85)
$\mathbf{R}=\mathbf{D}+\sum_{i=1}^{k}\mathbf{s}_{i}D^{-1}\cdot\mathbf{c}_{i}^{+},$

where $\mathbf{s}_{i}$ are symmetries, and $\mathbf{c}_{i}$ are cosymmetries
of the integrable system.

We will find vector quasi-local recursion and symplectic operator for
Boussinesq type systems from Theorem 2. For each of systems (61), (62) and
(64) we have found two Hamiltonian operators and present a bi-Hamiltonian for
of the system. For system (63) we also found two Hamiltonian operators.
However, we were able to write in the bi-Hamiltonian form a higher symmetry,
but not the system itself.

### System (62)

This system admits the symmetries $\mathbf{s}_{1}=\begin{bmatrix}[c]1\\\
0\end{bmatrix}$, $\mathbf{s}_{2}=\begin{bmatrix}[c]u_{1}\\\
v_{1}\end{bmatrix}$ and cosymmetries $\mathbf{c}_{1}=\begin{bmatrix}[c]0\\\
1\end{bmatrix}$, $\mathbf{c}_{2}=\begin{bmatrix}[r]-v_{1}\\\
u_{1}\end{bmatrix}$, variational derivatives of conserved densities $v$ and
$u_{1}v$ respectively. Using these symmetries and cosymmetries we can write
recursion, symplectic and Hamiltonian operators
$\mathbf{R},\mathbf{S}_{1},\mathbf{S}_{2},\mathbf{H}_{1}$ and $\mathbf{H}_{2}$
in quasi-local form as follows:

$\mathbf{R}=\begin{bmatrix}[c]\frac{3}{4}v&D\\\\[10.00002pt]
R_{21}&\frac{3}{4}v\end{bmatrix}+\frac{1}{4}\begin{bmatrix}u_{1}\\\
v_{1}\end{bmatrix}D^{-1}\cdot\begin{bmatrix}0&1\end{bmatrix}+\frac{1}{4}\begin{bmatrix}[c]1\\\
0\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]-v_{1}&u_{1}\end{bmatrix},$

where
$R_{21}=D^{5}+\frac{5}{2}u_{1}D^{3}+\frac{15}{4}u_{2}D^{2}+\left(\frac{9}{4}u_{3}+u_{1}^{2}\right)+\frac{1}{2}u_{4}+u_{1}u_{2}$,

$\mathbf{S}_{1}=\begin{bmatrix}[r]0&-1\\\
1&0\end{bmatrix},\qquad\mathbf{H}_{1}=\mathbf{S}_{1}^{-1}=\begin{bmatrix}[r]0&1\\\
-1&0\end{bmatrix}$

and

$\displaystyle\mathbf{S}_{2}=\mathbf{S}_{1}\mathbf{R}=\begin{bmatrix}[c]-R_{21}&-\frac{3}{4}v\\\\[15.00002pt]
\frac{3}{4}v&D\end{bmatrix}+\frac{1}{4}\begin{bmatrix}[r]-v_{1}\\\
u_{1}\end{bmatrix}D^{-1}\cdot\begin{bmatrix}0&1\end{bmatrix}+\frac{1}{4}\begin{bmatrix}[c]0\\\
1\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]-v_{1}&u_{1}\end{bmatrix},$
$\displaystyle\mathbf{H}_{2}=-\mathbf{R}\mathbf{H}_{1}=\begin{bmatrix}D&-\frac{3}{4}v\\\\[10.00002pt]
\frac{3}{4}v&-R_{21}\end{bmatrix}+\frac{1}{4}\begin{bmatrix}[c]u_{1}\\\
v_{1}\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]1&0\end{bmatrix}+\frac{1}{4}\begin{bmatrix}[c]1\\\
0\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]u_{1}&v_{1}\end{bmatrix}.$

The Hamiltonian operators $\mathbf{H}_{1}$ and $\mathbf{H}_{2}$ are compatible
[3]. A bi-Hamiltonian pencil can be obtained from $\mathbf{H}_{2}$ by shifting
$v\to v+\text{const}$ and two corresponding Hamiltonian forms of the system
are given by

$\begin{bmatrix}u\\\
v\end{bmatrix}_{t}=\mathbf{H}_{1}\frac{\delta}{\delta\mathbf{u}}\left(\frac{1}{2}v^{2}-\frac{1}{2}u_{2}^{2}+\frac{1}{3}u_{1}^{3}\right)=\mathbf{H}_{2}\frac{\delta}{\delta\mathbf{u}}\left(-2v\right).$

### System (64)

This system admits the symmetry $(1,0)$ and the conserved density
$\rho=v-u_{1}^{2}$, with variational derivative $(2u_{2},1)$, allowing us to
write the quasi-local expression for a recursion operator

$\mathbf{R}=\begin{bmatrix}[c]-u_{1}&0\\\
D^{3}+\left(2v-6u_{1}^{2}\right)D+v_{1}&3u_{1}\end{bmatrix}+\begin{bmatrix}[c]1\\\
0\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]2u_{2}&1\end{bmatrix}.$

Two symplectic operators are

(86) $\mathbf{S}_{1}=\begin{bmatrix}4u_{1}D+2u_{2}&-1\\\ 1&0\end{bmatrix}$

and

(87)
$\mathbf{S}_{2}=\mathbf{S}_{1}\mathbf{R}=\begin{bmatrix}-D^{3}+2(u_{1}^{2}-v)D+\left(2u_{1}u_{2}-v_{1}\right)&u_{1}\\\
-u_{1}&0\end{bmatrix}+\begin{bmatrix}2u_{2}\\\
1\end{bmatrix}D^{-1}\cdot\begin{bmatrix}2u_{2}&1\end{bmatrix}.$

Two Hamiltonian operators have the form

$\mathbf{H}_{1}=\mathbf{S}_{1}^{-1}=\begin{bmatrix}[c]\phantom{-}0&1\\\
-1&4u_{1}D+2u_{2}\end{bmatrix}$

and

$\mathbf{H}_{2}=-\mathbf{R}\mathbf{H}_{1}=\begin{bmatrix}[c]0&-3u_{1}\\\
3u_{1}&-D^{3}-\left(2v+6u_{1}^{2}\right)D-v_{1}-6u_{1}u_{2}\end{bmatrix}+\begin{bmatrix}[c]1\\\
0\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]1&0\end{bmatrix}$

These Hamiltonian operators are related by the argument shift $u\to
u+\text{const}\,x$. We can write the system in bi-Hamiltonian form as

$\begin{bmatrix}u\\\
v\end{bmatrix}_{t}=\mathbf{H}_{1}\frac{\delta}{\delta\mathbf{u}}\left(\frac{1}{2}v^{2}-\frac{1}{2}u_{2}^{2}-\frac{1}{2}u_{1}^{4}\right)=\mathbf{H}_{2}\frac{\delta}{\delta\mathbf{u}}\left(-u_{1}v+u_{\text{1}}^{3}\right).$

### System (63)

The simplest symmetries and cosymmetries of this system are

$\mathbf{s}_{1}=\begin{bmatrix}[c]1\\\ 0\end{bmatrix},\
\mathbf{s}_{2}=\begin{bmatrix}[c]v\\\ v_{t}\end{bmatrix},\
\mathbf{s}_{3}=\begin{bmatrix}[c]u_{1}\\\
v_{1}\end{bmatrix},\quad\mathbf{c}_{1}=\begin{bmatrix}[c]-2u_{2}\\\
1\end{bmatrix},\ \mathbf{c}_{2}=\begin{bmatrix}[c]-4u_{1}u_{2}-2v_{1}\\\
2u_{1}\end{bmatrix},$

where the latter two are the variational derivatives of conserved densities
$v+u_{1}^{2}$, and $\frac{2}{3}u_{1}^{3}+2u_{1}v$. Using those we find the
quasi-local recursion operator

$\mathbf{R}=\begin{bmatrix}[c]\frac{1}{2}u_{2}D+\frac{1}{2}u_{1}^{3}+\frac{3}{2}u_{1}v&D\\\\[10.00002pt]
R_{21}&\frac{5}{2}u_{2}D+\frac{1}{2}u_{1}^{3}+\frac{3}{2}u_{1}v+2u_{3}\end{bmatrix}\\\
{}+\frac{1}{2}\begin{bmatrix}v\\\
u_{4}+2u_{2}v+2u_{1}^{2}u_{2}\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]-2u_{2}&1\end{bmatrix}\\\
{}+\frac{1}{2}\begin{bmatrix}[c]u_{1}\\\
v_{1}\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]-2u_{1}u_{2}-v_{1}&2u_{1}\end{bmatrix}.$

with

$R_{21}=D^{5}+\frac{5}{2}(u_{1}^{2}+v)D^{3}+(\frac{15}{2}u_{1}u_{2}+2v_{1})D^{2}\\\
\textstyle{}+(u_{4}+2u_{1}^{2}v+v^{2}+\frac{9}{2}u_{1}u_{3}+4u_{2}^{2}+\frac{1}{2}v_{2})D\\\
\textstyle{}+\frac{1}{2}u_{1}^{2}v_{1}+2u_{1}u_{2}v+\frac{1}{2}vv_{1}+2u_{1}^{3}u_{2}+u_{1}u_{4}$

and the symplectic operator

$\mathbf{S}=\begin{bmatrix}[c]S_{11}&-D^{2}-\frac{3}{2}u_{1}^{2}-\frac{1}{2}v\\\
D^{2}+\frac{3}{2}u_{1}^{2}+\frac{1}{2}v&0\end{bmatrix}\\\
{}+\frac{1}{2}\begin{bmatrix}[c]u_{1}\\\
v_{1}\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]-2u_{2}&1\end{bmatrix}+\frac{1}{2}\begin{bmatrix}[c]1\\\
0\end{bmatrix}D^{-1}\cdot\begin{bmatrix}[c]-2u_{1}u_{2}-2v_{1}&2u_{1}\end{bmatrix}$

with

$S_{11}=-u_{1}D^{3}-\frac{3}{2}u_{2}D^{2}-\left(2u_{1}^{3}+2u_{1}v+\frac{1}{2}u_{3}\right)D-3u_{1}^{2}u_{2}-u_{2}v-u_{1}v_{1}.$

Two Hamiltonian operators are

$\mathbf{H}_{1}=\begin{bmatrix}0&u_{1}\\\
-u_{1}&-D^{3}-2(u_{1}^{2}+v)D-2u_{1}u_{2}-v_{1}\end{bmatrix}+\begin{bmatrix}[c]1\\\
0\end{bmatrix}D^{-1}\begin{bmatrix}[c]1&0\end{bmatrix},$

and

$\mathbf{H}_{2}=\begin{bmatrix}-2u_{1}D-u_{2}&H_{12}\\\
H_{21}&H_{22}\end{bmatrix}+\begin{bmatrix}[c]u_{1}\\\
v_{1}\end{bmatrix}D^{-1}\begin{bmatrix}[c]v&v_{t}\end{bmatrix}+\begin{bmatrix}[c]v\\\
v_{t}\end{bmatrix}D^{-1}\begin{bmatrix}[c]u_{1}&v_{1}\end{bmatrix}$

with

$\displaystyle\begin{split}H_{12}&=-2D^{4}-5(u_{1}^{2}+v)D^{2}-2(5u_{1}u_{2}+3v_{1})D\\\
&\quad{}-u_{1}^{4}-vu_{1}^{2}-5u_{3}u_{1}-2v^{2}-3u_{2}^{2}-2v_{2},\\\
H_{21}&=-H_{12}-2v_{1}D-v_{2},\\\
H_{22}&=2u_{1}D^{5}+5u_{2}D^{4}+(4u_{1}^{3}+2vu_{1}+16u_{3})D^{3}\\\
&\quad{}+(18u_{2}u_{1}^{2}+3v_{1}u_{1}+3vu_{2}+19u_{4})D^{2}\\\
&\quad{}+(-4vu_{1}^{3}+16u_{3}u_{1}^{2}-4v^{2}u_{1}+8u_{2}^{2}u_{1}+v_{2}u_{1}-6u_{2}v_{1}+7vu_{3}+10u_{5})D\\\
&\quad{}-2v_{1}u_{1}^{3}-6vu_{2}u_{1}^{2}+5u_{4}u_{1}^{2}-4vv_{1}u_{1}+6u_{2}u_{3}u_{1}-2u_{2}^{3}-2v^{2}u_{2}-v_{1}u_{3}\\\
&\quad{}-4u_{2}v_{2}+3vu_{4}+2u_{6}.\end{split}$

Using the first Hamiltonian operator we can write the system in the
Hamiltonian form:

$\begin{bmatrix}u\\\
v\end{bmatrix}_{t}=\mathbf{H}_{1}\frac{\delta}{\delta\mathbf{u}}\left(-u_{1}v-\frac{1}{3}u_{1}^{3}\right).$

We cannot find any conserved density to express the system with the operator
$\mathbf{H}_{2}$. However, possibly all higher symmetries of (63) can be
written in bi-Hamiltonian form. In particular, for the symmetry

$\mathbf{s}=\begin{bmatrix}[c]v_{2}+\frac{1}{2}u_{2}^{2}+u_{1}^{2}v+\frac{1}{2}v^{2}+\frac{1}{6}u_{1}^{4}\\\\[11.38109pt]
\begin{gathered}\textstyle
u_{6}+3u_{2}v_{2}+4u_{3}v_{1}+3vu_{4}+3u_{1}^{2}u_{4}+2u_{1}vv_{1}+\frac{2}{3}u_{1}^{3}v_{1}\\\
{}+12u_{1}u_{2}u_{3}+4u_{1}^{2}u_{2}v+2u_{2}v^{2}+2u_{1}^{4}u_{2}+4u_{2}^{3}\end{gathered}\end{bmatrix}$

we have

$\mathbf{s}=\mathbf{H}_{2}\frac{\delta\
}{\delta\mathbf{u}}\left(-\frac{v}{2}-\frac{u_{1}^{2}}{2}\right)=\mathbf{H}_{1}\frac{\delta\
}{\delta\mathbf{u}}\left(u_{2}v_{1}+\frac{1}{2}u_{1}u_{2}^{2}-\frac{2}{3}vu_{1}^{3}-\frac{1}{2}v^{2}u_{1}-\frac{u_{1}^{5}}{6}\right).$

###  8.1. System (61).

In this case, writing the operators in quasi-local form requires symmetries
and cosymmetries that explicitly depend on $x$ and $t$555When we expand the
operators to pseudo-differential series this dependence disappears.. Using the
symmetries

$\mathbf{s}_{1}=\begin{bmatrix}[c]1\\\ 0\end{bmatrix},\
\mathbf{s}_{2}=\begin{bmatrix}[c]x\\\ 0\end{bmatrix},\
\mathbf{s}_{3}=\begin{bmatrix}[c]t\\\ 1\end{bmatrix},\
\mathbf{s}_{4}=\begin{bmatrix}[c]xt\\\ x\end{bmatrix}$

and cosymmetries

$\mathbf{c}_{1}=\begin{bmatrix}[r]v_{3}\\\ -u_{3}\end{bmatrix},\
\mathbf{c}_{2}=\begin{bmatrix}[c]2u_{6}+4u_{2}u_{4}+4u_{3}^{2}\\\
-2v_{2}\end{bmatrix},\
\mathbf{c}_{3}=\begin{bmatrix}[c]-2tu_{3}^{2}-2v_{2}-2tu_{2}u_{4}-\frac{1}{2}xv_{3}-tu_{6}\\\
tv_{2}+\frac{1}{2}xu_{3}+u_{2}\end{bmatrix},$

the first two being variational derivatives of the conserved densities
$\rho_{1}=u_{2}v_{1}$ and $\rho_{2}=v_{1}^{2}-u_{3}^{2}+\frac{2}{3}u_{2}^{3}$,
we find the following recursion operator:

$\mathbf{R}=\begin{bmatrix}\frac{3}{4}v_{1}&D\\\\[5.69054pt]
D^{5}+\tfrac{5}{2}u_{2}D^{3}+\tfrac{5}{4}u_{3}D^{2}+\left(u_{4}+u_{2}^{2}\right)D-\tfrac{1}{2}u_{5}-u_{2}u_{3}&\frac{3}{4}v_{1}\end{bmatrix}\\\
+\frac{1}{4}\begin{bmatrix}[c]x\\\
0\end{bmatrix}D^{-1}\cdot\mathbf{c}_{1}^{+}+\frac{1}{4}\begin{bmatrix}[c]t\\\
1\end{bmatrix}D^{-1}\cdot\mathbf{c}_{2}^{+}+\frac{1}{2}\begin{bmatrix}[c]1\\\
0\end{bmatrix}D^{-1}\cdot\mathbf{c}_{3}^{+}.$

Two differential symplectic operators are

$\mathbf{S}_{1}=\begin{bmatrix}[c]0&-D^{2}\\\
D^{2}&0\end{bmatrix},\quad\mathbf{S}_{2}=\begin{bmatrix}[c]-D^{7}-\frac{5}{2}u_{2}D^{5}-\frac{25}{4}u_{3}D^{4}&-\frac{3}{4}v_{1}D^{2}-v_{2}D-\frac{1}{4}v_{3}\\\\[5.0pt]
\frac{3}{4}v_{1}D^{2}+\frac{1}{2}v_{2}D&D^{3}+\frac{1}{2}u_{2}D+\frac{1}{4}u_{3}\end{bmatrix}.$

The recursion operator found is
$\mathbf{R}=\mathbf{S}_{1}^{-1}\mathbf{S}_{2}$.

A Hamiltonian operator is

(88) $\mathbf{H}_{1}=\mathbf{S}_{1}^{-1}=\begin{bmatrix}[c]0&D^{-2}\\\
-D^{-2}&0\end{bmatrix}\\\ =\begin{bmatrix}x\\\
0\end{bmatrix}D^{-1}\begin{bmatrix}t&1\end{bmatrix}+\begin{bmatrix}t\\\
1\end{bmatrix}D^{-1}\begin{bmatrix}x&0\end{bmatrix}-\begin{bmatrix}1\\\
0\end{bmatrix}D^{-1}\begin{bmatrix}xt&x\end{bmatrix}-\begin{bmatrix}xt\\\
x\end{bmatrix}D^{-1}\begin{bmatrix}1&0\end{bmatrix}$

and the system can be written in Hamiltonian form as

$\begin{bmatrix}[c]u\\\
v\end{bmatrix}_{t}=\frac{1}{2}\mathbf{H}_{1}\frac{\delta\rho_{2}}{\delta{\bf
u}}.$

Another Hamiltonian operator is

$\mathbf{H}_{2}=\mathbf{R}\mathbf{H}_{1}=\begin{bmatrix}0&0\\\
0&-D^{3}-\frac{5}{2}u_{2}D-\frac{5}{4}u_{3}\end{bmatrix}-\frac{1}{4}\begin{bmatrix}x\\\
0\end{bmatrix}D^{-1}\begin{bmatrix}u_{1}&v_{1}\end{bmatrix}-\frac{1}{4}\begin{bmatrix}u_{1}\\\
v_{1}\end{bmatrix}D^{-1}\begin{bmatrix}x&0\end{bmatrix}\\\
{}-\frac{1}{2}\begin{bmatrix}t\\\
1\end{bmatrix}D^{-1}\begin{bmatrix}v&v_{t}\end{bmatrix}-\frac{1}{2}\begin{bmatrix}v\\\
v_{t}\end{bmatrix}D^{-1}\begin{bmatrix}t&1\end{bmatrix}+\frac{1}{4}\begin{bmatrix}1\\\
0\end{bmatrix}D^{-1}\begin{bmatrix}xu_{1}+2tv&xv_{1}+2v+2tv_{t}\end{bmatrix}\\\
{}+\frac{1}{4}\begin{bmatrix}1+xu_{1}+2tv\\\
xv_{1}+2v+2tv_{t}\end{bmatrix}D^{-1}\begin{bmatrix}1&0\end{bmatrix}.$

The given Hamiltonian operators are compatible since they are related by the
shift $v\to v+\text{const}\,x$. Nevertheless, $\mathbf{H}_{2}$ cannot be used
to write a second Hamiltonian form of the system, in a local form. Observe
that, using the substitution $u\to Du$, $v\to Dv$ transforms (61) into system
(62).

##  9\. Concluding remarks

We have explained a derivation of explicit integrability expressions for
certain degenerate systems, that we are defined as regularly diagonalisable
system. These systems are those that can be formally diagonalized after the
usual diagonalisation procedure for non-degenerate systems. We have used the
explicit integrability conditions for classifying a subclass of Boussinesq-
type integrable systems, showing that non-trivial systems arise even in the
simplest case.

##  10\. Acknowledgements

The authors want to thank V.S. Novikov for insightful comments. RHH was
supported by projects PID2021-124473NB-I00 and PID2019-106802GB-I00 from the
Ministry of Science and Innovation of Spain, and is a member of the
International Society of Nonlinear Mathematical Physics is$\sim$mp. VVS was
supported by the FAPESR grant 2023/14060-7 (Brazil). He is grateful to IME USP
for its hospitality during his visit to Sao Paulo. Both authors want to
express their gratitude to the Institut des Hautes Études Scientifiques IHES
(France) for their hospitality and invaluable support.

## References

  * [1] Mikhailov, A V, Shabat, A B and Yamilov, R I, Symmetry approach to classification of nonlinear equations. Complete lists of integrable systems. Russian Math. Surveys, 42(4), 1–63, 1987.
  * [2] Mikhailov A V, Shabat A B and Sokolov V V, Symmetry Approach to Classification of Integrable Equations, What is integrability? Ed. V.E. Zakharov, Springer Series in Nonlinear Dynamics, Springer-Verlag, 115–184, 1991.
  * [3] Olver P J, Applications of Lie Groups to Differential Equations, (2nd edn), Graduate Texts in Mathematics, 107, Springer-Verlag, New York, 1993.
  * [4] Zhiber A V and Shabat A B, Klein-Gordon equations with a nontrivial group, Sov. Phys. Dokl., 24(8), 608–609, 1979.
  * [5] Fokas A S, Symmetries and integrability, Studies in Applied Mathematics, 7̱7(3), 253–299, 1987.
  * [6] Abellanas L and Galindo A, Evolution equations with high order conservation laws. J. Math. Phys., 24, 3, 504–509, 1983.
  * [7] Mikhailov A V, Novikov V S and Wang J P. On Classification of Integrable Nonevolutionary Equations. Studies in Applied Mathematics, 118(4), 419–457, 2007.
  * [8] Novikov V S and Wang J P, Symmetry structure of integrable nonevolutionary equations. Studies in Applied Mathematics, 119(4), 393–428, 2007.
  * [9] Hernández Heredero R, Sokolov V V, and Shabat A B, A new class of linearizable equations. J. Phys. A: Math. Gen. 36, L605–L614, 2003.
  * [10] Ibragimov N Kh and Shabat A B, Infinite Lie-Bäklund algebras, Funct. anal. appl. 14(4), 313–315, 1980.
  * [11] Sokolov V V and Shabat A B, Classification of integrable evolution equations, Soviet Scientific Reviews, Section C, 4, 221–280, 1984.
  * [12] Chen H H, Lee Y C and Liu C S, Integrability of nonlinear Hamiltonian systems by inverse scattering method, Phys. Scr., 20(3-4), 490–492, 1979.
  * [13] Svinolupov S I, Sokolov V V, Evolution equations with nontrivial conservative laws, Funct. anal. appl., 16(4), 317–319, 1982.
  * [14] Sokolov V V, Hamiltonian property of the Krichever-Novikov equation, Sov. Math. Dokl., 30, 44–47, 1984.
  * [15] Maltsev A Ya and Novikov S P, On the local systems Hamiltonian in the weakly non-local Poisson brackets, Phys. D, 156(1-2), 53–80, 2001.
  * [16] Demskoy D K and Sokolov V V, On recursion operators for elliptic models, Nonlinearity, 21, 1253–1264, 2008.
  * [17] Sanders J A and Wang J P, Integrable systems and their recursion operators. Nonlinear Analysis: Theory, Methods & Applications, 47(8), 5213–5240, 2001.
  * [18] Caparrós Quintero A and Hernández Heredero R, Formal recursion operators of integrable nonevolutionary equations and Lagrangian Systems. J. Phys. A: Math. Theor. 51, 385201, 2018.
  * [19] Dorfman I Ya, Dirac Structures and Integrability of Nonlinear Evolution Equations, John Wiley&Sons, Chichester, 1993.