On Recursion Operators for Full-Fledged Nonlocal Symmetries of the Reduced Quasi-classical Self-dual Yang–Mills Equation

Abstract

We introduce the idea of constructing recursion operators for full-fledged nonlocal symmetries and apply it to the reduced quasi-classical self-dual Yang–Mills equation. It turns out that the discovered recursion operators can be interpreted as infinite-dimensional matrices of differential functions which act on the generating vector functions of the nonlocal symmetries simply by matrix multiplication. To the best of our knowledge, there are no other examples of such recursion operators in the literature so far, so our approach is completely innovative. Further, we investigate the algebraic properties of the discovered operators and discuss the ${\mathbb {R}}$-algebra structure on the set of all recursion operators for full-fledged nonlocal symmetries of the equation in question. Finally, we illustrate the action of the obtained recursion operators on particularly chosen full-fledged symmetries and emphasize their advantages compared to the action of traditionally used recursion operators for shadows.

1 Introduction

Partial differential equations play an important role in physics and other applied sciences because they describe many natural phenomena. Therefore, it is no wonder that they have been intensively investigated over the past few decades. A very important tool in the study of systems of nonlinear partial differential equations is the symmetry analysis: having the symmetry of the system and one of its solutions, we can construct another one. A Lie subalgebra $\mathfrak g$ of the Lie algebra of all infinitesimal symmetries of the given system can be used to determine special types of solutions—the ${\mathfrak {g}}$-invariant solutions; for many nonlinear systems, these solutions are the only exact solutions that we are able to obtain (see, e.g., [18, 28]); thus, they are very important from the mathematical as well as physical point of view. Systems that possess infinitely many generalized symmetries are often predisposed to linearization. This can be typically achieved either through a change of variables or by employing the inverse scattering method, see, e.g., [1, 10, 28, 33].

A significant role in the study of symmetries of systems of partial differential equations is played by recursion operators, because they allow us to generate infinite hierarchies of symmetries. Roughly speaking, given a differential equation, if we know one of its symmetries and if we have a recursion operator at our disposal, we can apply this operator to generate a new symmetry; then we apply the recursion operator to this new symmetry and so on, arriving at the infinite hierarchy of symmetries of the equation in question.

In the context of symmetry transformations, recursion operators were first introduced by P. J. Olver in 1977 (see [27]) and further used and developed by other authors, see, e.g., [5, 14, 28] and references therein. Within this approach, the recursion operators arise as pseudodifferential operators, i.e., differential operators expressed as the linear combinations of total derivative operators and their formal inverses, the coefficients being the (matrix valued) differential functions. Hence, in general, they may contain the expressions like $D_x^{-1}$, and sometimes it may happen that it is not clear how these formal operators should act on a given symmetry.

An alternative view on the recursion operators, based on the concept of Bäcklund transformations, was proposed in [29], and this idea was further developed and given the rigorous mathematical background by Marvan in [24]. The general idea of this geometric approach is that a recursion operator can be viewed as a Bäcklund auto-transformation for the linearized equation, i.e., the equation defining the (generating functions, or characteristics of) symmetries of the original nonlinear PDE.

Treating a recursion operator as a Bäcklund auto-transformation, it can be applied to any symmetry, but the result may depend on new additional quantities. More precisely, for a given partial differential system ${\mathcal {E}}$, it is natural to consider its differential covering $\tilde{{\mathcal {E}}}$, which can be roughly viewed as an overdetermined system that introduces additional dependent variables (nonlocal variables, pseudopotentials) such that ${\mathcal {E}}$ implies the compatibility conditions of $\tilde{{\mathcal {E}}}$, and subsequently to study the recursion operators which act within the set of $\tilde{{\mathcal {E}}}$-shadows of nonlocal symmetries of ${\mathcal {E}}$. The recursion operators constructed in this way may provide us with infinite hierarchies of $\tilde{{\mathcal {E}}}$-shadows. Nevertheless, concerning full-fledged nonlocal symmetries the things are not so simple (see, e.g., [6]): not every $\tilde{{\mathcal {E}}}$-shadow can be lifted to a full-fledged nonlocal symmetry of ${\mathcal {E}}$ (i.e., a symmetry of $\tilde{{\mathcal {E}}}$) in the covering $\mathcal {{\tilde{E}}}$. Hence, the recursion operators computed in this way (and usually they are computed in this way) do not provide us immediately with the infinite-dimensional Lie algebra of symmetries. However, the infinite-dimensional Lie algebras of nonlocal symmetries of ${\mathcal {E}}$ play an important role in the theory of integrable systems and can be used for the study of the latter (see, e.g., [17]). So, it would be very useful to know a technique that would enable us to introduce a sufficiently rich differential covering and subsequently construct recursion operators that map (in this covering) full-fledged nonlocal symmetries to other ones.

One of the often studied properties of recursion operators is their Nijenhuis torsion (see, e.g., [17] and references therein for more details). Recall that if a recursion operator ${\mathcal {R}}$ is Nijenhuis, i.e., its Nijenhuis torsion $N_{{\mathcal {R}}}$ is trivial, and $\varphi _1$ and $\varphi _2$ are two commuting symmetries such that the Lie derivatives $L_{\varphi _1}({\mathcal {R}})$ and $L_{\varphi _2}({\mathcal {R}})$ vanish, then the symmetries $\mathcal R^k(\varphi _1)$ and ${\mathcal {R}}^l(\varphi _2)$ mutually commute for all $k,l \in {\mathbb {N}}$. Let us point out that it is still not completely clear (see, e.g., [18]) how one should verify the Nijenhuis property for general recursion operators acting on nonlocal symmetries of PDEs in higher dimensions. However, the Nijenhuis torsion of the recursion operators is defined in terms of the Jacobi brackets of unspecified symmetries and their images under the action of the operator in question. To compute these Jacobi brackets, it is not enough to work with the shadows only but we are forced to know full-fledged symmetries. Hence, it seems that the full-fledged recursion operators proposed in this paper could also be a good stepping stone for investigation in this direction.

In the present paper, we develop the idea of constructing recursion operators for full-fledged nonlocal symmetries and demonstrate it on the reduced quasi-classical self-dual Yang–Mills equation (rYME)

$$\begin{aligned} u_{yz} = u_{tx} - u_z u_{xx}+u_xu_{xz}. \end{aligned}$$

(1)

To the best of our knowledge, this equation was first introduced in [12] and subsequently studied, e.g., in [15, 21, 25]; it also appears in the context of the classification of integrable four-dimensional systems associated with sixfolds in Gr(4,6) as the equation arising as the compatibility condition for a canonical form of general linearly degenerate systems, see [9]. Our approach to the problem of construction of full-fledged recursion operators consists in exploiting certain dependencies inherent in the “common” recursion operators that act only on the shadows of nonlocal symmetries. Based on certain observations, we present an explicit construction of the recursion operator as a standard mapping on the set of full-fledged (in the given covering) nonlocal symmetries of the rYME (1). However, let us point out that this our demonstration example was chosen more or less randomly and in principle we could select another one. Indeed, based on specific computational experiments performed on other equations, we are pretty sure that the stated explicit dependency between recursion operators below can also be discovered for other 4D linearly degenerate Lax-integrable equations, namely for all equations listed in Tab. 2 in [9] (which include the rYME, as well as, for example, the 4D Martínez-Alonso Shabat equation, see, e.g., [32] and references therein for more details), or for the modified Martínez-Alonso Shabat equation [26] (cf. also [2]). Our approach is surely also applicable to 3D equations like the (modified) Veronese web equation, the rdDym equation, the universal hierarchy equation and the Pavlov equation, see, e.g., [3, 11, 16, 22, 34] for more information about these equations, although in this case the set of possible recursion operators seems to be (quite logically) much poorer (and hence less interesting) then in the case of 4D equations.

The structure of this paper is as follows: in Sect. 1, we recall a few necessary notions and facts from the geometrical theory of PDEs, for example the notions of linearization, symmetry, differential covering, tangent equation and recursion operator.

In Sect. 2, we proceed to the construction of the differential covering and a full-fledged recursion operator for symmetries of the rYME (1). In order to construct a sufficiently large covering, we use techniques that have been used successfully before, e.g., in [3, 16, 20] or [32]. They are based on the existence of the Lax representation of the equation in question which contains two non-removable parameters $\lambda$ and $\mu$: we expand the nonlocal variable into a formal Laurent series in $\lambda$ and in $\mu$, respectively, and arrive at three inequivalent infinite-dimensional coverings. Their Whitney product forms the desired covering which enables us to find an explicit dependence between the components of full-fledged nonlocal symmetries whose shadows are related by a ‘common’ recursion operator. The resulting full-fledged recursion operator for the rYME (1) has surprisingly nice properties: its action on the generating vector functions of nonlocal symmetries can be viewed just as (infinite-dimensional) matrix multiplication, the corresponding matrix having only finitely many nonzero entries in each row and its entries being only nonlocal differential functions living in the given Whitney product.

Section 3 is devoted to the investigation of algebraic properties of the above constructed recursion operator. It turns out that it is an automorphism of the ${\mathcal {A}}$-module ${\mathcal {A}}^{{\mathbb {N}}}$ (where ${\mathcal {A}}$ is a suitably chosen ring of nonlocal functions living in the Whitney product) whose inverse is also a full-fledged recursion operator. This implies that the presented recursion operator is an isomorphism on the real vector space of full-fledged nonlocal symmetries whose generating functions lie in ${\mathcal {A}}^{{\mathbb {N}}}$. We also find another invertible full-fledged recursion operator that possesses the same nice properties as the first one. We prove that all the just obtained recursion operators and their inverses pair-wise mutually commute, which implies that the set of all the compositions of these recursion operators carries a commutative group structure with respect to the composition operation. Finally, we conjecture that the two just found full-fledged recursion operators together with their inverses are the generators of the whole associative algebra of all recursion operators for the rYME (1).

In Sect. 4, we highlight the advantages of the full-fledged recursion operators for symmetries of the rYME (1) compared to the traditionally used recursion operators for shadows and discuss other benefits and possibilities associated with a recursion operator for full-fledged symmetries. To this end, we demonstrate with one example how a recursion operator for shadows acts on a particular (shadow of) symmetry and describe the challenges faced in constructing a hierarchy of full-fledged symmetries using such operators. Subsequently, we present several examples demonstrating how easily (in contrast to the previous) the application of the discovered ‘full-fledged’ recursion operators to particular full-fledged symmetries can be carried out, and we highlight that in this way one can quite straightforwardly obtain the full-fledged forms of those nonlocal symmetries which would be otherwise only with great difficulty constructed from their shadows.

2 Preliminaries

In this section, we briefly summarize some basic notions and facts from nonlocal geometry of PDEs which are crucial for the next sections of this paper. For more detailed exposition, we refer the reader to the classical book [6], the overview paper [17] and/or, e.g., [18, 19].

Let $\pi :{\mathbb {R}}^m\times {\mathbb {R}}^n\rightarrow {\mathbb {R}}^n$ be a trivial bundle where $x^1,\dots ,x^n$ are coordinates standing for the independent variables in ${\mathbb {R}}^n$ and the coordinates $u^1,\dots ,u^m$ in ${\mathbb {R}}^m$ denote the dependent variables. Let $\pi _{\infty }:J^{\infty }(\pi )\rightarrow {\mathbb {R}}^n$ be the bundle of infinite jets, the coordinates in $J^{\infty }(\pi )$ being $x^i,u^j,u_{\sigma }^j$, where $i=1,\dots ,n$, $j=1,\dots ,m$, and $\sigma$ are symmetric multi-indices of arbitrary finite length $|\sigma |$ that consist of the integers $1,\dots ,n$; the coordinate $u^j_\sigma$ stands for the derivative of the dependent variable $u^j$ with respect to the independent variables stated in the multi-index $\sigma$ here. The infinite jet space $J^{\infty }(\pi )$ is equipped with the total derivatives operators

$$\begin{aligned} D_{x^i}\equiv D_i = \frac{\partial }{\partial x^i} + \sum \limits _{j,\sigma } u_{\sigma i}^j \frac{\partial }{\partial u_\sigma ^j},\quad i=1,\dots ,n, \end{aligned}$$

that provide us with the so-called Cartan distribution ${\mathcal {C}}$ on $J^{\infty }(\pi )$. Below $D_{\sigma }$ denotes the composition of the total derivative operators with respect to those independent variables that are stated in the multi-index $\sigma$.

A differential function on $J^{\infty }(\pi )$ is a smooth function F on $J^{\infty }(\pi )$ of finite order

$$\begin{aligned} \textrm{ord}\ F=\textrm{max}\left\{ |\sigma |\ |\ \partial F/\partial u^j_\sigma \ne 0\right\} . \end{aligned}$$

The ${\mathbb {R}}$-algebra of differential functions is denoted by ${\mathcal {F}}(\pi )$. Each function $F\in {\mathcal {F}}(\pi )$ of order k can be identified with a nonlinear differential operator $\Delta _F(f^1,\dots ,f^m)=F(x,f^j,\ldots , \frac{\partial ^{|\sigma |}f^j}{\partial x^{\sigma }}, \ldots )$ of order k; thus, each system of differential functions $\{ F^{\alpha }\} \subset {\mathcal {F}}(\pi )$ provides us with a system of differential equations

$$\begin{aligned} \Delta _{F^{\alpha }}(f^1,\dots ,f^m)=0, \end{aligned}$$

where $f=(f^1(x),\dots ,f^j(x),\dots ,f^m(x))$ are the unknown functions in the independent variables $x=(x^1,\ldots , x^n)$, and vice versa.

Let $F=(F^1,\dots ,F^{r})$ be a vector differential function on $J^{\infty }(\pi )$, i.e., $F^{\alpha }\in {\mathcal {F}}(\pi )$ for each $\alpha =1,\dots ,r$. From now on, an (infinitely prolonged) differential equation ${\mathcal {E}}$ is the submanifold

$$\begin{aligned} {\mathcal {E}}\equiv {\mathcal {E}}_F=\left\{ \theta \in J^{\infty }(\pi )\ |\ D_{\sigma }(F^{\alpha })(\theta )=0\ \forall \alpha ,\sigma \right\} \subset J^{\infty }(\pi ), \end{aligned}$$

(2)

(it corresponds to the system of differential equations given by the functions $F^{\alpha }$ considered together with all its differential consequences).The structure of the equation ${\mathcal {E}}$ is given by the Cartan distribution restricted to ${\mathcal {E}}$ which is spanned by the $\pi _{\infty }$-horizontal vector fields

$$\begin{aligned} {{\bar{D}}}_i={D_i}\big |_{{\mathcal {E}}},\quad i=1,\dots ,n. \end{aligned}$$

A solution of ${\mathcal {E}}$ is a section $s=(u^1=s^1(x),\dots ,u^m=s^m(x))$ of the bundle $\pi$ such that its infinite jet $j_{\infty }(s)=(u^1_{\sigma }=\frac{\partial ^{|\sigma |}s^1}{\partial x^{\sigma }},\dots ,u^m_{\sigma }=\frac{\partial ^{|\sigma |}s^1}{\partial x^{\sigma }})$ (here $\sigma$ runs over all multi-indices) lies in ${\mathcal {E}}$. The solutions correspond to the maximal n-dimensional integral manifolds of the Cartan distribution ${\mathcal {C}}$.

Suppose that it is possible to choose and fix a set ${\mathbb {I}}$ of internal coordinates on ${\mathcal {E}}$ (this can be done by locally resolving the equations $F^1=F^2=\dots =F^r=0$ with respect to certain partial derivatives). Denote ${{\mathcal {F}}}({\mathcal {E}})$ the algebra of smooth functions on $J^{\infty }(\pi )$ restricted to ${\mathcal {E}}$ and expressed in the internal coordinates. A $\pi _{\infty }$-vertical vector field

$$\begin{aligned} S=\sum _{u_\sigma ^j\in {\mathbb {I}}}S_{\sigma }^j\frac{\partial }{\partial u^j_{\sigma }}, \quad S_{\sigma }^j\in {{\mathcal {F}}}({\mathcal {E}}), \end{aligned}$$

is a (higher infinitesimal) symmetry of the equation ${\mathcal {E}}$ if it holds $[{\mathcal {C}},S]\subset {\mathcal {C}}$. It can be shown that symmetries of the equation ${\mathcal {E}}$ are exactly the evolutionary vector fields

$$\begin{aligned} {\textbf{E}}_{\Phi } = \sum \limits _{u_\sigma ^j \in {\mathbb {I}}} \bar{D}_{\sigma }(\varphi _j)\frac{\partial }{\partial u_\sigma ^j}, \quad \varphi _j\in {\mathcal {F}}({\mathcal {E}}), \end{aligned}$$

(3)

where $\Phi =(\varphi _1,\varphi _2, \ldots , \varphi _m)$ satisfies the condition

$$\begin{aligned} \ell _{{\mathcal {E}}} (\Phi ) \equiv {\ell _F}\big |_{{\mathcal {E}}} (\Phi ) = 0, \end{aligned}$$

(4)

and

$$\begin{aligned} \ell _F=\begin{pmatrix} \sum _\sigma \frac{\partial F^1}{\partial u_\sigma ^1}D_\sigma &{} \sum _\sigma \frac{\partial F^1}{\partial u_\sigma ^2}D_\sigma &{} \ldots \\[3mm] \sum _\sigma \frac{\partial F^2}{\partial u_\sigma ^1}D_\sigma &{} \sum _\sigma \frac{\partial F^2}{\partial u_\sigma ^2}D_\sigma &{} \ldots \\ \vdots &{} \vdots &{} \ddots \end{pmatrix} \end{aligned}$$

denotes the operator of linearization of F. The vector-valued differential function $\Phi \in \left( {\mathcal {F}}({\mathcal {E}})\right) ^m$ is called the generating function of the symmetry in question.

The set of all symmetries of the equation ${\mathcal {E}}$ equipped with the commutator forms a Lie algebra ${{\,\textrm{sym}\,}}({\mathcal {E}})$. It is isomorphic to the Lie algebra of all solutions to (4), the Lie algebra structure being given by the Jacobi bracket $\left\{ \cdot ,\cdot \right\}$ defined by the formula:

$$\begin{aligned} \left\{ \Phi _1,\Phi _2 \right\} = {\textbf{E}}_{\Phi _1}(\Phi _2)-\textbf{E}_{\Phi _2}(\Phi _1). \end{aligned}$$

(5)

Due to this isomorphism, we will not distinguish below between the symmetries and their generating functions.

A recursion operator for symmetries of ${\mathcal {E}}$ (in the classical sense of Olver, see, e.g., [28]) is an ${\mathbb {R}}$- linear operator ${\mathcal {R}}:({\mathcal {F}}({\mathcal {E}}))^m\rightarrow ({\mathcal {F}}({\mathcal {E}}))^m$ with the property that whenever $\Phi$ is a generating function of a symmetry of ${\mathcal {E}}$, so is ${\mathcal {R}}(\Phi )$.

Consider the space $\tilde{{\mathcal {E}}} = {\mathbb {R}}^l \times \mathcal E$, where $1\le l \le \infty$ and ${\mathcal {E}}$ is determined by (2). Denote the coordinates in ${\mathbb {R}}^l$ as $w^\eta$ and assume that the space $\tilde{{\mathcal {E}}}$ is endowed with the operators

$$\begin{aligned} {{\tilde{D}}}_i = {{\bar{D}}}_i + \sum \limits _\eta X_i^\eta \mathchoice{\frac{\partial }{\partial w^\eta }}{\partial /\partial w^\eta }{\partial /\partial w^\eta }{\partial /\partial w^\eta },\quad X_i^{\eta }\in {\mathcal {F}}(\tilde{{\mathcal {E}}}), \end{aligned}$$

(6)

which pair-wise commute, i.e.,

$$\begin{aligned} {[}{{\tilde{D}}}_i, {{\tilde{D}}}_j]={{\bar{D}}}_i(X_j^\eta )-{{\bar{D}}}_j(X_i^\eta ) + \sum \limits _\nu \left( X_i^\nu \mathchoice{\frac{\partial X_j^\eta }{\partial w^\nu }}{\partial X_j^\eta /\partial w^\nu }{\partial X_j^\eta /\partial w^\nu }{\partial X_j^\eta /\partial w^\nu } - X_j^\nu \mathchoice{\frac{\partial X_i^\eta }{\partial w^\nu }}{\partial X_i^\eta /\partial w^\nu }{\partial X_i^\eta /\partial w^\nu }{\partial X_i^\eta /\partial w^\nu } \right) = 0 \end{aligned}$$

(7)

for all i, j and $\eta$ (in other words, the Cartan distribution $\tilde{{\mathcal {C}}}$ on $\tilde{{\mathcal {E}}}$ given by these vector fields is integrable). We say that the projection $\tau : \tilde{{\mathcal {E}}} \rightarrow {\mathcal {E}}$ to the second factor of $\tilde{{\mathcal {E}}}$ is a (differential) covering over the equation ${\mathcal {E}}$ and the coordinates $w^\eta$ are called nonlocal coordinates or nonlocal variables in this covering. The covering space ${\tilde{{\mathcal {E}}}}$ can be equivalently regarded as a differential equation $\tilde{{\mathcal {E}}}={\mathcal {E}}\cap {\mathcal {E}}_{\tau }\subset J^{\infty }({\tilde{\pi }})$, where the equation ${\mathcal {E}}_{\tau }$ is given by the relations

$$\begin{aligned} {\mathcal {E}}_{\tau }:\ \frac{\partial w^{\eta }}{\partial x^i}=X_i^\eta (x,u_{\sigma }^j,w^{\gamma }), \quad i=1,\ldots ,n, \; \eta =1,\ldots ,l, \end{aligned}$$

(8)

where $w^1,\dots ,w^l$ are considered to be new unknown functions of the variable $x=(x^1,\dots ,x^n)$. The equation ${\mathcal {E}}_{\tau }$ is called the covering equation or the Lax representation for ${\mathcal {E}}$. In the case when the system (8) is given in the form of two equations for each nonlocal variable, it is referred to as a Lax pair. The compatibility conditions for the system (8) are differential consequences of the equation ${\mathcal {E}}$.

Symmetries of the equation $\tilde{{\mathcal {E}}}$ are said to be nonlocal symmetries of the equation ${\mathcal {E}}$ in this covering. Thus, (full-fledged) nonlocal $\tau$-symmetries of ${\mathcal {E}}$ are ${{{\tilde{\pi }}}}^{\infty }$-vertical vector fields X on $\tilde{{\mathcal {E}}}$ that preserve the Cartan distribution $\tilde{{\mathcal {C}}}$ on $\tilde{{\mathcal {E}}}$, i.e., $[X,\tilde{{\mathcal {C}}}]\subset \tilde{{\mathcal {C}}}$, which means that X must be exactly of the form

$$\begin{aligned} X\equiv {\textbf{E}}_{(\Phi ,\Psi )}=\sum \limits _{u_\sigma ^j\in {\mathbb {I}}} \tilde{D}_{\sigma }(\varphi _j)\frac{\partial }{\partial u_\sigma ^j} +\sum \limits _\eta \psi ^\eta \mathchoice{\frac{\partial }{\partial w^\eta }}{\partial /\partial w^\eta }{\partial /\partial w^\eta }{\partial /\partial w^\eta }, \end{aligned}$$

(9)

where $\Phi =(\varphi _1, \varphi _2,\ldots ,\varphi _m)$, $\Psi =(\psi _1, \psi _2, \ldots ,\psi _l)$ are vector-valued functions on $\tilde{{\mathcal {E}}}$ satisfying the system

$$\begin{aligned} {{\tilde{\ell }}}_{{\mathcal {E}}}(\Phi )&=0, \end{aligned}$$

(10)

$$\begin{aligned} {{\tilde{D}}}_i(\psi ^\eta )&= {{\tilde{\ell }}}_{X_i^\eta }(\Phi ) + \sum \limits _\nu \mathchoice{\frac{\partial X_i^\eta }{\partial w^\nu }}{\partial X_i^\eta /\partial w^\nu }{\partial X_i^\eta /\partial w^\nu }{\partial X_i^\eta /\partial w^\nu }\psi ^\nu , \end{aligned}$$

(11)

and ${{\tilde{\ell }}}_{{\mathcal {E}}}$ denotes the natural lift of the differential operator ${\ell }_{{\mathcal {E}}}$ from ${\mathcal {E}}$ to $\tilde{{\mathcal {E}}}$. Solutions of the Eq. (10) are called shadows of nonlocal symmetries. Let us point out that some authors call the shadows of nonlocal symmetries shortly as nonlocal symmetries. However, whenever we talk about $\tau$-nonlocal symmetry in the following, we always mean a full-fledged $\tau$-nonlocal symmetry.

Let $\tau _1:\tilde{{\mathcal {E}}}_1\rightarrow {\mathcal {E}}$ and $\tau _2:\tilde{{\mathcal {E}}}_2\rightarrow {\mathcal {E}}$ be two coverings over ${\mathcal {E}}$. Denote $v^1,v^2,\ldots ,v^k$ the nonlocal variables in $\tau _1$ and $w^1,w^2,\ldots ,w^l$ the nonlocal variables in $\tau _2$, $1\le k,l \le \infty$. Let the Cartan distributions $\tilde{{\mathcal {C}}}^{1}$ on $\tilde{{\mathcal {E}}}_1$, resp. $\tilde{{\mathcal {C}}}^{2}$ on $\tilde{{\mathcal {E}}}_2$ be given by the vector fields ${{\tilde{D}}}_i^{(1)}={{\bar{D}}}_i+X_i^{(1)}$, resp. $\tilde{D}_i^{(2)}={{\bar{D}}}_i+X_i^{(2)}$, $i=1,\dots ,n$, where $X_i^{(1)}= \sum _{\eta }{V_i^{\rho }}\frac{\partial }{\partial v^{\rho }}$ and $X_i^{(2)}= \sum _{\eta }{W_i^{\eta }}\frac{\partial }{\partial w^{\eta }}.$ Consider the direct product $\tilde{{\mathcal {E}}}_1\times \tilde{{\mathcal {E}}}_2$ and its subset $\tilde{{\mathcal {E}}}_1\oplus \tilde{{\mathcal {E}}}_2=\left\{ (y_1,y_2)\in \tilde{{\mathcal {E}}}_1\times \tilde{{\mathcal {E}}}_2 \ |\ \tau _1(y_1)=\tau _2(y_2)\right\}$ equipped with the (integrable) distribution ${{\mathcal {C}}}^{\oplus }$ defined by the fields $\tilde{D}_i^{\oplus }={{\tilde{D}}}_i+ X_i^{(1)}+X_i^{(2)}$. Then the projection $\tau _1\oplus \tau _2:\tilde{{\mathcal {E}}}_1\oplus \tilde{{\mathcal {E}}}_2\rightarrow {\mathcal {E}}, (\tau _1\oplus \tau _2)(y_1,y_2)=\tau _1(y_1)=\tau _2(y_2)$ is called the Whitney product of $\tau _1$ and $\tau _2$. The corresponding covering equation is given by the relations

$$\begin{aligned} {\mathcal {E}}_{\tau _1\oplus \tau _2}:\frac{\partial v^{\rho }}{\partial x^i}&=V_i^\rho (x,u_{\sigma }^j,v^{\gamma _1}),\\ \frac{\partial w^{\eta }}{\partial x^i}&=W_i^\eta (x,u_{\sigma }^j,w^{\gamma _2}),\!\!\!&i&=1,\ldots ,n, \\ {}{} & {} \rho ,\gamma _1&=1,\dots ,k, \quad \eta ,\gamma _2=1,\ldots ,l. \end{aligned}$$

The tangent equation $\mathcal{T}\mathcal{E}\subset J^{\infty }({\bar{\pi }})$ of the equation ${\mathcal {E}}$ is given by the system

$$\begin{aligned} \mathcal{T}\mathcal{E}:\ \begin{array}{r} F^\alpha \left( x,\ldots , u^j_{\sigma }, \ldots \right) =0,\\ \ell _F(p) = 0, \end{array} \end{aligned}$$

where $p=(p^1,p^2,\ldots ,p^m)$ are new dependent variables, the projection $t:{\mathcal {T}}{\mathcal {E}}\rightarrow {\mathcal {E}}, (x,u_{\sigma }^j,p_{\sigma })\mapsto (x,u_{\sigma })$ is the so-called tangent covering of ${\mathcal {E}}$. The sections of ${\mathcal {T}}{\mathcal {E}}$ that map the Cartan distribution on ${\mathcal {E}}$ to the Cartan distribution on ${\mathcal {T}}{\mathcal {E}}$ are the symmetries of ${\mathcal {E}}$.

A Bäcklund transformation between the equations ${\mathcal {E}}_1$ and ${\mathcal {E}}_2$ is the diagram

where $\tau _1$ and $\tau _2$ are coverings. If ${\mathcal {E}}_1={\mathcal {E}}_2={\mathcal {E}}$, then the Bäcklund transformation of ${\mathcal {E}}$ is called the Bäcklund auto-transformation. Any Bäcklund auto-transformation of $\mathcal{T}\mathcal{E}$ relates (shadows of) symmetries of ${\mathcal {E}}$ to each other, which means that it can be interpreted as a recursion operator. We discuss this approach widely in the next Sect. 2, where we describe our idea of constructing of recursion operators for full-fledged nonlocal symmetries.

3 Coverings and the Recursion Operator for Full-Fledged Nonlocal Symmetries of the rYME

Consider the rYME (1), where x, y, z, t are the independent variables and u stands for the dependent variable. The manifold ${\mathcal {E}}$ corresponding to this equation lies in the jet space $J^{\infty }(\pi )$, where $\pi :{\mathbb {R}}\times {\mathbb {R}}^4\rightarrow {\mathbb {R}}^4$. Then the internal coordinates on ${\mathcal {E}}$ are $u_{x^it^j},u_{x^iy^kt^j},u_{x^iz^lt^j}$ where $i,j\ge 0,\ k,l>0$, and the Cartan distribution on ${\mathcal {E}}$ is spanned by the total derivative operators $D_x,D_y,D_z,D_t$ on $J^{\infty }(\pi )$ restricted to ${\mathcal {E}}$, that is

As it has already been mentioned in the introductory part of this paper, the rYME (1) was introduced for the first time in [12], including its one-parameter isospectral Lax representation

$$\begin{aligned} s_t=u_zs_x+\varkappa s_z, \quad s_y=(u_x+\varkappa )s_x, \end{aligned}$$

(12)

cf. also [9]. However, it turns out that the Lax pair (12) can be easily modified by the substitution $\varkappa \rightarrow -\frac{z+\mu }{t+\lambda }$ into the Lax pair which contains two parameters and explicitly depends on the independent variables t, z. Indeed, it is straightforward to verify that the equations

$$\begin{aligned} w_t=u_zw_x-\frac{z+\mu }{t+\lambda }\,w_z, \quad w_y=u_xw_x-\frac{z+\mu }{t+\lambda }\,w_x, \end{aligned}$$

(13)

where $\lambda ,\mu \in {\mathbb {R}}$ are spectral parameters and w is a nonlocal variable, also define the Lax representation of the rYME (1).

Remark 1

In fact, the “nonlocal variables” produced by Lax pairs for equations in more than two independent variables are not nonlocal variables in the sense of the definition stated in Sect. 1. However, as we will discuss below, these “nonlocal variables” provide us with an infinite number of nonlocal variables for the given equation. Thus, until we analyze this problem for the case of the rYME (1) in more details, we will call w and other variables given by the Lax pairs for ${\mathcal {E}}$ to be just nonlocal variables.

To obtain coverings of (1), we consider the following scheme applied before, for example, in [3, 16, 32], cf. also [23, 30]:

• (1) The nonlocal variable w given by the Lax pair (13) is expanded in a formal series in the spectral parameter $\lambda$, resp. $\mu$:

  $$\begin{aligned} w=\sum \limits _{i=-\infty }^\infty \lambda ^i w_i, \; \mathrm {resp.}\quad w=\sum \limits _{i=-\infty }^\infty \mu ^i w_i. \end{aligned}$$

• (2) By substituting this expansion into the Lax pair (13) and eliminating the parameter $\lambda$, resp. $\mu$, we get an infinite series of nonlocal quantities $w_i$, $i \in {\mathbb {Z}}$.

• (3) To achieve the proper definition of nonlocal variables $w_i$, we consider two reductions of the series $w_i$: (a) the negative reduction $w_i=0$ for $i>0$ and (b) the positive reduction $w_i=0$ for $i<0$.

• (4) In this way, we arrive at the following three inequivalent infinite-dimensional coverings $\tau ^q$, $\tau ^m$ and $\tau ^r$ (given in the form of Lax pairs) with the nonlocal variables $q_\alpha ,m_\alpha ,r_\beta$, respectively, where

  $$\begin{aligned} \tau ^q:\tilde{{\mathcal {E}}}^q\rightarrow {\mathcal {E}}&\quad \left| \begin{array}{l} q_{-1} = 0,\\ q_{\alpha ,t} = u_zq_{\alpha ,x} + q_{\alpha -1,z},\\[2pt] q_{\alpha ,y} = u_xq_{\alpha ,x} + q_{\alpha -1,x}, \quad \alpha \ge 0; \end{array}\right. \end{aligned}$$

  (14)
  $$\begin{aligned} \tau ^m:\tilde{{\mathcal {E}}}^m\rightarrow {\mathcal {E}}&\quad \left| \begin{array}{l} m_{-1} = 0,\\ m_{\alpha ,t} = u_z m_{\alpha ,x}-\displaystyle \frac{z}{t}m_{\alpha ,z}-\frac{1}{t}m_{\alpha -1,z},\\[3mm] m_{\alpha ,y} = u_x m_{\alpha ,x}-\displaystyle \frac{z}{t}m_{\alpha ,x}-\frac{1}{t}m_{\alpha -1,x}, \quad \alpha \ge 0; \end{array}\right. \end{aligned}$$

  (15)

  and

  $$\begin{aligned} \tau ^r:\tilde{{\mathcal {E}}}^r\rightarrow {\mathcal {E}}&\quad \left| \begin{array}{l} r_{-1} = x,\quad r_0=-u\\ r_{\beta ,x} = u_xr_{\beta -1,x} - r_{\beta -1,y},\\[2pt] r_{\beta ,z} = u_zr_{\beta -1,x} - r_{\beta -1,t}, \quad \beta \ge 1. \end{array}\right. \end{aligned}$$

  (16)

Remark 2

Let us note that the coverings $\tau ^q$ and $\tau ^m$ are non-Abelian, whereas $\tau ^r$ is Abelian. Thus, the nonlocal variables $r_\beta$, $\beta \ge 1$, are interconnected with infinite series of two-component nonlocal conservation laws of the rYME (1).

Remark 3

Whereas the coverings $\tau ^q$ and $\tau ^r$ can also be derived from the one-parameter Lax pair (12), the covering $\tau ^m$ can be obtained only from the two-parameter Lax representation (13).

As noted in Remark 1, Eqs. (14), (15) and (16) are not in the form of the covering Eq. (8). In order to identify the nonlocal coordinates provided by these formulas, we have to introduce the infinite number of formal variables $q_{\alpha }^{a,b}$, $m_{\alpha }^{a,b}$ and $r_{\beta }^{a,b},$ $a,b,=0,1,\dots$, where $q_{\alpha }^{0,0}=q_{\alpha }$, $m_{\alpha }^{0,0}=m_{\alpha }$, $r_{\beta }^{0,0}=r_{\beta }$ and

$$\begin{aligned} q_{\alpha ,x}^{a,b}= & {} q_{\alpha }^{a+1,b}, \quad q_{\alpha ,z}^{a,b}=q_{\alpha }^{a,b+1},\quad q_{\alpha ,y}^{a,b}=(u_xq_{\alpha }^{1,0} + q_{\alpha -1}^{1,0})_{x^{a}z^b},\\ q_{\alpha ,t}^{a,b}= & {} (u_zq_{\alpha }^{1,0} + q_{\alpha -1}^{0,1})_{x^{a}z^b},\\ m_{\alpha ,x}^{a,b}= & {} m_{\alpha }^{a+1,b},\quad m_{\alpha ,z}^{a,b}=m_{\alpha }^{a,b+1},\quad m_{\alpha ,y}^{a,b}=\left( u_x m_{\alpha }^{1,0}-\displaystyle \frac{z}{t}m_{\alpha }^{1,0}-\frac{1}{t}m_{\alpha -1}^{1,0}\right) _{x^{a}z^b},\\ m_{\alpha ,t}^{a,b}= & {} \left( u_z m_{\alpha }^{1,0}-\displaystyle \frac{z}{t}m_{\alpha }^{0,1}-\frac{1}{t}m_{\alpha -1}^{0,1}\right) _{x^{a}z^b},\\ r_{\beta ,y}^{a,b}= & {} r_{\beta }^{a+1,b},\quad r_{\beta ,t}^{a,b}=r_{\beta }^{a,b+1},\\ r_{\beta ,x}^{a,b}= & {} (u_xr_{\beta -1,x} - r_{\beta -1}^{1,0})_{y^at^b},\quad r_{\beta ,z}^{a,b}=( u_zr_{\beta -1,x} - r_{\beta -1}^{0,1})_{y^at^b}. \end{aligned}$$

Note that all the right-hand sides of the equalities above are (even though sometimes inductively) expressible in terms of the variables $q_{\alpha }^{a,b}$, $m_{\alpha }^{a,b}$ and $r_{\beta }^{a,b},$ thus, they provide us with the covering equations ${\mathcal {E}}_{\tau ^q}$, ${\mathcal {E}}_{\tau ^m}$ and ${\mathcal {E}}_{\tau ^r}$ in the form (8), $q_{\alpha }^{a,b}$, $m_{\alpha }^{a,b}$ and $r_{\beta }^{a,b},$ $\alpha =0,1,\dots ,\beta =1,2,\dots$, being the true nonlocal variables.

The total derivative operators in the covering $\tau ^q$ now take the form (cf. (6))

$$\begin{aligned} {\tilde{D}}_x^{(q)}={{\bar{D}}}_x\!+\!X^{(q)},\!\quad {\tilde{D}}_y^{(q)}={{\bar{D}}}_y\!+\!Y^{(q)},\!\quad {\tilde{D}}_z^{(q)}={{\bar{D}}}_z\!+\!Z^{(q)},\!\quad {\tilde{D}}_t^{(q)}={{\bar{D}}}_t+T^{(q)}, \end{aligned}$$

where

$$\begin{aligned} X^{(q)}&=\sum _{\alpha =1}^{\infty }\sum _{a,b=0}^{\infty }q_{\alpha }^{a+1,b}\frac{\partial }{\partial q_{\alpha }^{a,b}},\qquad&Y^{(q)}=\sum _{\alpha =1}^{\infty }\sum _{a,b=0}^{\infty }(u_xq_{\alpha ,x} + q_{\alpha -1,x})_{x^{a}z^b}\frac{\partial }{\partial q_{\alpha }^{a,b}},\\ Z^{(q)}&=\sum _{\alpha =1}^{\infty }\sum _{a,b=0}^{\infty }q_{\alpha }^{a,b+1}\frac{\partial }{\partial q_{\alpha }^{a,b}},&T^{(q)}=\sum _{\alpha =1}^{\infty }\sum _{a,b=0}^{\infty }(u_zq_{\alpha ,x} + q_{\alpha -1,z})_{x^{a}z^b}\frac{\partial }{\partial q_{\alpha }^{a,b}}. \end{aligned}$$

The total derivative operators $D_x^{(m)}$, $D_x^{(r)}$, $D_y^{(m)}$, $D_y^{(r)}$, etc. in the coverings $\tau ^m$ and $\tau ^r$ are given in a similarly way.

In what follows, we will consider the Whitney product $\tau ^W=\tau ^q \oplus \tau ^m \oplus \tau ^r$ of all three coverings and carry out all calculations in $\tau ^W$. The covering space ${\tilde{\mathcal {E}}}^W={\mathcal {E}}_{\tau ^q \oplus \tau ^m \oplus \tau ^r}\cap {\mathcal {E}}$ can be considered to be an equation in the jet space $J^{\infty }({{\tilde{\pi }}})$ given by the rYME (1) and Eqs. (14)–(16) for $\alpha \ge 0$, $\beta \ge 1$, the coordinates in $J^{\infty }({{\tilde{\pi }}})$ being $x,y,z,t,u_{\sigma },q_{\alpha ,\sigma }^{a,b},m_{\alpha ,\sigma }^{a,b},r_{\beta ,\sigma }^{a,b},\ a,b\ge 0$. The Cartan distribution on $\tilde{{\mathcal {E}}}^W$ is spanned by the total derivative operators

$$\begin{aligned} {{\tilde{D}}}_x&={{\bar{D}}}_x+X^{(q)}+X^{(m)}+X^{(r)},&{{\tilde{D}}}_y&=\bar{D}_y+Y^{(q)}+Y^{(m)}+Y^{(r)},\\ {{\tilde{D}}}_z&={{\bar{D}}}_z+Z^{(q)}+Z^{(m)}+Z^{(r)},&{{\tilde{D}}}_t&=\bar{D}_t+T^{(q)}+T^{(m)}+T^{(r)}. \end{aligned}$$

As follows directly from (9)–(11), every full-fledged nonlocal symmetry of ${\mathcal {E}}$ in the Whitney product $\tau ^W$ (for shortness we will write only $\tau ^W$-symmetry below) is of the form

$$\begin{aligned}{} & {} {\textbf{E}}_P ={\textbf{E}}^W_{p_0}+\sum _{\alpha =0}^{\infty }\sum _{a,b=0}^{\infty } \\{} & {} \quad \left( {{\tilde{D}}}_x^a{{\tilde{D}}}_z^b(p_{\alpha }^{q})\frac{\partial }{\partial q_{\alpha }^{a,b}} +{{\tilde{D}}}_x^a{{\tilde{D}}}_z^b(p_{\alpha }^{m})\frac{\partial }{\partial m_{\alpha }^{a,b}}+{{\tilde{D}}}_y^a{{\tilde{D}}}_t^b(p_{\alpha +1}^{r})\frac{\partial }{\partial r_{\alpha +1}^{a,b}}\right) , \end{aligned}$$

where

and $p_0, p_{\alpha }^{q},p_{\alpha }^{m}, p_1^{r}, p_{\beta }^{r}$, $\alpha \ge 0,\beta \ge 2$, are smooth functions on $\tilde{{\mathcal {E}}}^W$ satisfying the conditions:

$$\begin{aligned} {{\tilde{D}}}_{yz}(p_0) - {{\tilde{D}}}_{tx}(p_0) + u_z {{\tilde{D}}}_{xx}(p_0)-u_x{{\tilde{D}}}_{xz}(p_0)+u_{xx}{{\tilde{D}}}_z(p_0)-u_{xz}{{\tilde{D}}}_x(p_0)&=0, \end{aligned}$$

(NS1)

$$\begin{aligned} {{\tilde{D}}}_t(p_\alpha ^q) - u_z{{\tilde{D}}}_x(p_\alpha ^q) - q_{\alpha ,x}{{\tilde{D}}}_z(p_0)-{{\tilde{D}}}_z(p_{\alpha -1}^q)&=0, \end{aligned}$$

(NS2)

$$\begin{aligned} {{\tilde{D}}}_y(p_\alpha ^q) -u_x{{\tilde{D}}}_x(p_\alpha ^q) - q_{\alpha ,x}{{\tilde{D}}}_x(p_0)-{{\tilde{D}}}_x(p_{\alpha -1}^q)&=0, \end{aligned}$$

(NS3)

$$\begin{aligned} {{\tilde{D}}}_t(p_\alpha ^m) -u_z{{\tilde{D}}}_x(p_\alpha ^m) - m_{\alpha ,x}{{\tilde{D}}}_z(p_0)+\frac{z}{t}{{\tilde{D}}}_z(p_\alpha ^m)+\frac{1}{t}{{\tilde{D}}}_z(p_{\alpha -1}^m)&=0, \end{aligned}$$

(NS4)

$$\begin{aligned} {{\tilde{D}}}_y(p_\alpha ^m) -u_x{{\tilde{D}}}_x(p_\alpha ^m) - m_{\alpha ,x}{{\tilde{D}}}_x(p_0)+\frac{z}{t}{{\tilde{D}}}_x(p_\alpha ^m)+\frac{1}{t}{{\tilde{D}}}_x(p_{\alpha -1}^m)&=0, \end{aligned}$$

(NS5)

$$\begin{aligned} {{\tilde{D}}}_x(p_1^r) - {{\tilde{D}}}_y(p_0)+2u_x{{\tilde{D}}}_x(p_0)&=0, \end{aligned}$$

(NS6)

$$\begin{aligned} {{\tilde{D}}}_z(p_1^r) -{{\tilde{D}}}_t(p_0)+u_x{{\tilde{D}}}_z(p_0)+u_z{{\tilde{D}}}_x(p_0)&=0, \end{aligned}$$

(NS7)

$$\begin{aligned} {{\tilde{D}}}_x(p_\beta ^r) - u_x{{\tilde{D}}}_x(p_{\beta -1}^r)-r_{\beta -1,x}{{\tilde{D}}}_x(p_0)+D_y(p_{\beta -1}^r)&=0, \end{aligned}$$

(NS8)

$$\begin{aligned} {{\tilde{D}}}_z(p_\beta ^r) - u_z{{\tilde{D}}}_x(p_{\beta -1}^r)-r_{\beta -1,x}{{\tilde{D}}}_z(p_0)+D_t(p_{\beta -1}^r)&=0, \end{aligned}$$

(NS9)

where $p_{-1}^q=p_{-1}^m=0$.

Thus, any $\tau ^W$-symmetry can be (uniquely) represented by the vector-valued generating function $P=\left[ p_0, p_0^q, p_0^m, p_1^r,p_1^q, p_1^m, p_2^r,p_2^q, p_2^m\ldots \right]$. From now on, unless otherwise stated, by a nonlocal symmetry of ${\mathcal {E}}$ we mean a generating function of a $\tau ^W$-symmetry of ${\mathcal {E}}$, and $\textrm{sym}^{\tau ^W}({\mathcal {E}})$ will denote the set of all generating functions of $\tau ^W$-symmetries of ${\mathcal {E}}$ (equipped with the Lie algebra structure given by the Jacobi bracket).

Remark 4

Let us note that we have used the formal nonlocal variables $q_{\alpha }^{a,b},m_{\alpha }^{a,b},r_{\beta }^{a,b}$ in order to justify the general form of the (generating function of) $\tau ^W$-symmetry of ${\mathcal {E}}$. However, since the equalities $q_{\alpha ,x^az^b}\equiv q_{\alpha ,x^az^b}^{0,0}=q_{\alpha }^{a,b},m_{\alpha ,x^az^b}\equiv m_{\alpha ,x^az^b}^{0,0}=m_{\alpha }^{a,b},r_{\beta ,y^at^b}\equiv r_{\beta ,y^at^b}^{0,0}=r_{\beta }^{a,b}$ hold on $\tilde{{\mathcal {E}}}^{W}$, we will use the more common notation $q_{\alpha ,x^az^b},m_{\alpha ,x^az^b},r_{\beta ,y^at^b}$ for them. Let us also emphasize that we sometimes use terms like $m_{\alpha ,t}, r_{\beta ,x}$, etc., in the explicit formulas presented below, even though $m_{\alpha ,t}, r_{\beta ,x}$, etc., are not internal coordinate on $\tilde{{\mathcal {E}}}^W.$ In fact, such terms are used to rewrite the huge resulting formulas into a more concise form.

The task to find a recursion operator for $\tau ^W$-symmetries for ${\mathcal {E}}$ now coincides with the problem to find a Bäcklund auto-transformation of the tangent equation ${\mathcal {T}} \tilde{{\mathcal {E}}}^W$ given by Eqs. (1) and (NS1)–(NS9) with the new dependent variables $p_{0,\sigma },p_{\alpha ,x^a z^b}^q, p_{\alpha ,x^a z^b}^m, p_{\beta ,y^a t^b}^r$, the tangent covering being

$$\begin{aligned}&t^W: {\mathcal {T}} {\tilde{{\mathcal {E}}}}^W \rightarrow {\tilde{\mathcal {E}}}^W,\\&\quad (x,y,z,t,u_\sigma ^j,q_{\alpha ,x^a z^b},m_{\alpha ,x^a z^b}, r_{\beta ,y^a t^b}, p_{0,\sigma },p_{\alpha ,x^a z^b}^q, p_{\alpha ,x^a z^b}^m, p_{\beta ,y^a t^b}^r)\\&\qquad \mapsto (x,y,z,t,u_\sigma ^j,q_{\alpha ,x^a z^b},m_{\alpha ,x^a z^b},r_{\beta ,y^a t^b}). \end{aligned}$$

The most straightforward way is to look for a Bäcklund auto-transformation expressible in the terms of the existing variables, i.e., to consider the covering equation common for both the copies of ${\mathcal {T}}{\tilde{{\mathcal {E}}}}^W$ to be just ${\mathcal {T}}{\tilde{{\mathcal {E}}}}^W$ itself. Due to the ${\mathbb {R}}$-linearity requirement imposed on the sought operator, we look for the desired relationship between the solutions to ${\mathcal {T}}{\tilde{{\mathcal {E}}}}^W$ in the following form: the components of the ‘new’ solution to ${\mathcal {T}}{\tilde{{\mathcal {E}}}}^W$ are to be (finite) linear combinations of the components $p_{0,\sigma },p_{\alpha ,x^a z^b}^q, p_{\alpha ,x^a z^b}^m, p_{\beta ,y^a t^b}^r$ of the ‘old’ solution to ${\mathcal {T}}{\tilde{{\mathcal {E}}}}^W$, the coefficients being arbitrary functions of the internal coordinates on ${\tilde{\mathcal {E}}}^W$. This ‘symmetry’, if it exists and we are able to describe it in an usable form, provides us with a recursion operator for $\tau ^W$-symmetries of ${\mathcal {E}}$.

The issue described above can be executed in the following way: we suppose that the quantity $P=\left[ p_0, p_0^q, p_0^m, p_1^r,p_1^q, p_1^m, p_2^r,p_2^q, p_2^m,\ldots \right]$ is a $\tau ^W$-symmetry of the rYME (1), i.e., the components of P satisfy the conditions (NS1)–(NS9) on ${\tilde{{\mathcal {E}}}}^W$, and subsequently try to solve Eqs. (NS1)–(NS9) with respect to the unknown quantity ${{\hat{P}}}=\left[ {{\hat{p}}}_0, {{\hat{p}}}_0^q, {{\hat{p}}}_0^m, {{\hat{p}}}_1^r,{{\hat{p}}}_1^q, {{\hat{p}}}_1^m, {{\hat{p}}}_2^r, {{\hat{p}}}_2^q, {{\hat{p}}}_2^m,\ldots \right]$ under the assumption that the components of ${{\hat{P}}}$ are linear combinations of the components of P and their derivatives up to some finite order, the coefficients being (possibly) arbitrary functions of the internal coordinates on $\tilde{{\mathcal {E}}}^W$.

Note that even though the idea of how to find a full-fledged recursion operator is quite straightforward, the implementation of the latter is not simple at all. The main difficulty of this construction consists in the fact that the recursion operator we are looking for is to act on an infinite sequence resulting in another infinite sequence, so that our success depends on whether we are able to find general formulas describing the relationships between all the components of the pre-image and its image.

One but not the only one result (see Sect. 3 below) we obtained within this approach is the following:

Proposition 1

Let $P=\left[ p_0, p_0^q, p_0^m, p_1^r,p_1^q, p_1^m, p_2^r,p_2^q, p_2^m,\ldots \right]$ be a $\tau ^W$-symmetry of the rYME (1), and let ${{\hat{P}}}=\left[ {{\hat{p}}}_0, {{\hat{p}}}_0^q, {{\hat{p}}}_0^m, {{\hat{p}}}_1^r,{{\hat{p}}}_1^q, {{\hat{p}}}_1^m, {{\hat{p}}}_2^r, {{\hat{p}}}_2^q, {{\hat{p}}}_2^m,\ldots \right]$ be a vector-valued function such that for each $\alpha \ge 0, \beta \ge 1$ it holds

$$\begin{aligned} {{\hat{p}}}_0&= u_x p_0+p_1^r, \end{aligned}$$

(17)

$$\begin{aligned} {{\hat{p}}}_\alpha ^q&= q_{\alpha ,x} p_0 + p_{\alpha -1}^q, \end{aligned}$$

(18)

$$\begin{aligned} {{\hat{p}}}_\alpha ^m&= m_{\alpha ,x} p_0 -\frac{1}{t} p_{\alpha -1}^m - \frac{z}{t} p_\alpha ^m , \end{aligned}$$

(19)

$$\begin{aligned} {{\hat{p}}}_\beta ^r&= r_{\beta ,x} p_0 - p_{\beta +1}^r. \end{aligned}$$

(20)

Then ${{\hat{P}}}$ is also a $\tau ^W$-symmetry of the rYME (1).

Proof

The proof is done by straightforward computations and is presented in Appendix A. $\square$

Remark 5

According to Proposition 1, Eqs. (17)–(20) provide us with a relation which maps any known $\tau ^W$-symmetry P to a new $\tau ^W$-symmetry ${{\hat{P}}}$. Hence, they can be viewed as the determining equations of a recursion operator for $\tau ^W$-symmetries of the rYME (1). In what follows, we will denote this recursion operator as ${\mathcal {R}}^q$. The notation will be clarified a bit later in Remark 10, see Sect. 3 below.

Remark 6

One can also easily observe that Eqs. (17)–(20) can be rewritten in the form

$$\begin{aligned} {{{\hat{P}}}}^T={\mathcal {R}}^q\cdot P^{T}, \end{aligned}$$

where ${}^T$ denotes the transpose of the row vectors P and ${{\hat{P}}}$, respectively, $\cdot$ is the usual matrix multiplication, and the recursion operator ${\mathcal {R}}^q$ is considered to be an infinite-dimensional matrix whose structure is illustrated by its $13 \times 16$ left-upper corner as follows:

$$\begin{aligned} {{\mathcal {R}}^q}= \left( \begin{array}{lllllllllllllllll} u_x &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ q_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ m_{0,x} &{} 0 &{} -z/t &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ r_{1,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} -1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ q_{1,x} &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}0 &{} 0 &{}\ldots \\ m_{1,x} &{} 0 &{} -1/t &{} 0 &{} 0 &{} -z/t &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ r_{2,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} -1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ q_{2,x} &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ m_{2,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} -1/t &{} 0 &{} 0 &{} -z/t &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ r_{3,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} -1 &{} 0 &{} 0 &{} 0 &{}\ldots \\ q_{3,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ m_{3,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} -1/t &{} 0 &{} 0 &{} -z/t &{} 0 &{} 0 &{} 0 &{} 0 &{}\ldots \\ r_{4,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} -1 &{}\ldots \\ \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \ddots \end{array}\right) \end{aligned}$$

As it can be easily verified, each column of the matrix ${\mathcal {R}}^q$ is itself a $\tau ^W$-symmetry of the rYME (1). However, this is perfectly natural, because if the matrix ${\mathcal {R}}^q$ is supposed to represent a recursion operator, then the images of the simplest nontrivial $\tau ^W$-symmetries $(1,0,0,\ldots ), (0,1,0,\ldots ),$ etc. must again be $\tau ^W$-symmetries.

Remark 7

Let us also note that using Eqs. (NS6)–(NS7) we can readily eliminate the quantity $p_1^r$ from Eq. (17). In this way, we obtain the conventional form (i.e., the form most frequently used in the current literature, cf., e.g., [2, 3, 16, 31, 32]) of the recursion operator ${\mathcal {R}}^q_{sh}$ for shadows of $\tau ^W$-symmetries of the rYME (1), which reads

$$\begin{aligned} \begin{aligned} {\mathcal {R}}^q_{sh}: {{\hat{p}}}_{0,x}&= p_{0,y}-u_x p_{0,x} + u_{xx} p_0,\\ {{\hat{p}}}_{0,z}&= p_{0,t}-u_z p_{0,x} + u_{xz} p_0. \end{aligned} \end{aligned}$$

(21)

As far as we know, even this recursion operator for shadows has not been, quite surprisingly, presented in literature yet.

4 Properties of the Recursion Operator ${\mathcal {R}}^q$ and Other Recursion Operators for $\tau ^W$-Symmetries of the rYME

In order to study the recursion operator ${\mathcal {R}}^q$ in more details, it is necessary to determine the ${\mathbb {R}}$-algebra ${\mathcal {A}}={\mathcal {A}}(\tilde{{\mathcal {E}}}^W)$ of functions (of the internal variables) on $\tilde{{\mathcal {E}}}^W$ that we take into consideration as the components of the generating functions of the $\tau ^W$-symmetries of $\tilde{{\mathcal {E}}}^W$.

According to the general theory (see, e.g., [6]), working in the coordinates in $J^{\infty }({{\tilde{\pi }}})$, the components of the generating functions of $\tau ^W$-symmetries are required to arise as the restrictions of differential functions on $J^{\infty }({{\tilde{\pi }}})$ to $\tilde{{\mathcal {E}}}^W$. However, even functions smooth on $J^{\infty }({\tilde{\pi }})$ expressed in the internal coordinates on $\tilde{{\mathcal {E}}}^W$ may not (as the functions of the internal coordinates) remain smooth globally any more. Thus, working in the internal coordinates on $\tilde{{\mathcal {E}}}^W$, it is reasonable to require the functions from ${\mathcal {A}}$ to be smooth only locally. For our purposes, it will be convenient to consider ${\mathcal {A}}$ to be the ${\mathbb {R}}$-algebra generated by the union ${\mathcal {L}}\cup {\mathcal {F}}$, where ${\mathcal {F}}$ denotes the algebra of smooth functions of the internal coordinates on $\tilde{{\mathcal {E}}}^W$ and ${\mathcal {L}}$ is the algebra of all Laurent polynomial functions¹ of the internal variables that depend on only finitely many of these variables (we will refer to the just defined functions as the differential Laurent polynomial functions on $\tilde{{\mathcal {E}}}^W$ below). Note that in case we are interested just in polynomial symmetries we can put ${\mathcal {L}}$ instead of ${\mathcal {A}}$ everywehere below and all the results will still remain valid.

Let ${\mathcal {A}}^{{\mathbb {N}}}=({{\mathcal {A}}(\tilde{{\mathcal {E}}}^W)})^{{\mathbb {N}}}$ denote the real vector space of all sequences of admissible functions on the equation $\tilde{{\mathcal {E}}}^W$, let ${\textrm{sym}}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})={\textrm{sym}}_{{\mathcal {A}}}(\tilde{{\mathcal {E}}}^W)\subset {\mathcal {A}}^{{\mathbb {N}}}$ denote the Lie algebra of the generating functions of full-fledged nonlocal $\tau ^W$-symmetries for (1) whose components lie in ${{{\mathcal {A}}}}$. The mapping ${\mathcal {R}}^q$ introduced in Sect. 2 can be readily seen to be a mapping ${\mathcal {A}}^{{\mathbb {N}}}\rightarrow {\mathcal {A}}^{{\mathbb {N}}}$ that is ${\mathbb {R}}$-linear and leaves the vector subspace $\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ of ${\mathcal {A}}^{{\mathbb {N}}}$ invariant, i.e., ${{{\mathcal {R}}}}^q(\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}}))\subset \textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$. Thus, it is a recursion operator in the sense of the standard definition presented in Sect. 1. Moreover, since ${{{\mathcal {R}}}}^q$ acts on ${\mathcal {A}}^{{\mathbb {N}}}$ just as the matrix multiplication operator described above, it can be easily concluded that ${{\mathcal {R}}}^q$ is also an ${\mathcal {A}}$-linear mapping:

Corollary 1

The recursion operator ${\mathcal {R}}^q$ is an ${\mathcal {A}}$-module ${{\mathcal {A}}}^{{\mathbb {N}}}$ endomorphism that leaves the set $\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ invariant.

It turns out that there exists another recursion operator for $\tau ^W$-symmetries of the rYME (1) within the class of ${\mathcal {A}}$-linear mappings.

Proposition 2

Let ${\mathcal {I}}$ be the identity mapping on ${\mathcal {A}}^{{\mathbb {N}}}$. Then the ${\mathcal {A}}$-linear mapping ${\mathcal {R}}^m:{\mathcal {A}}^{{\mathbb {N}}}\rightarrow {\mathcal {A}}^{{\mathbb {N}}}$ defined by the formula

$$\begin{aligned} {\mathcal {R}}^m = t\, {\mathcal {R}}^q + z\, {\mathcal {I}}, \end{aligned}$$

(22)

is a recursion operator for $\tau ^W$-symmetries of the rYME (1).

Proof

The proof is done by straightforward computations and is presented in Appendix B. $\square$

Remark 8

Note that the recursion operator ${\mathcal {R}}^m$ can be viewed as a matrix operator too, since the identity mapping ${\mathcal {I}}$ can be represented by the infinite-dimensional unit matrix. The $10 \times 14$ left-upper corner of the corresponding matrix is presented in Appendix C. Similarly, as in the case of ${\mathcal {R}}^q$ (see Remark 6) we can observe that each column of the matrix ${\mathcal {R}}^m$ is itself a $\tau ^W$-symmetry of the rYME (1).

Remark 9

As it can be seen immediately from the definition of the operator ${\mathcal {R}}^m$ (cf. (17) and (22)), the shadow $p_0$ of a symmetry $P \in {\textrm{sym}}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ and the shadow ${{\bar{p}}}_0$ of its ${\mathcal {R}}^m$-image ${{\bar{P}}} = {\mathcal {R}}^m(P)$ are related by the formula

$$\begin{aligned} {{\bar{p}}}_0 = (tu_x + z) p_0 + t p_1^r. \end{aligned}$$

(23)

As in the case of the recursion operator ${\mathcal {R}}^q$ (see Remark 7), we can eliminate the quantity $p_1^r$ from (23), which will provide us with the conventional form of the recursion operator ${\mathcal {R}}^m_{{\mathrm sh}}$ for shadows of $\tau ^W$-symmetries of the rYME (1):

$$\begin{aligned} \begin{aligned} {\mathcal {R}}^m_{sh}: {{\bar{p}}}_{0,x}&= t(p_{0,y}-u_x p_{0,x} + u_{xx} p_0)+zp_{0,x},\\ {{\bar{p}}}_{0,z}&= t(p_{0,t}-u_z p_{0,x} + u_{xz} p_0)+zp_{0,z}+p_0. \end{aligned} \end{aligned}$$

(24)

Proposition 3

The ${\mathcal {A}}$-module endomorphisms ${{{\mathcal {R}}}}^q:{{\mathcal {A}}}^{{\mathbb {N}}}\rightarrow {{\mathcal {A}}}^{{\mathbb {N}}}$, resp. ${{{\mathcal {R}}}}^m:{{\mathcal {A}}}^{{\mathbb {N}}}\rightarrow {{\mathcal {A}}}^{{\mathbb {N}}}$ are bijective, and their inverses $({{{\mathcal {R}}}}^{q})^{-1}$, resp. $({{{\mathcal {R}}}}^{m})^{-1}$ are given by the formulas

$$\begin{aligned} p_0&= \frac{{{\hat{p}}}_0^q}{q_{0,x}}, \end{aligned}$$

(25)

$$\begin{aligned} p_\alpha ^q&= -\frac{q_{\alpha +1,x}}{q_{0,x}} {{\hat{p}}}_0^q + {{\hat{p}}}_{\alpha +1}^q, \end{aligned}$$

(26)

$$\begin{aligned} p_\alpha ^m&= \sum \limits _{j=0}^\alpha \frac{t{{\hat{p}}}_j^m}{(-z)^{\alpha -j+1}}-\frac{t {{\hat{p}}}_0^q}{q_{0,x}} \sum \limits _{j=0}^\alpha \frac{m_{j,x}}{(-z)^{\alpha -j+1}}, \end{aligned}$$

(27)

$$\begin{aligned} p_\beta ^r&= \frac{r_{\beta -1,x}}{q_{0,x}}{{\hat{p}}}_0^q - {{\hat{p}}}_{\beta -1}^r, \end{aligned}$$

(28)

resp.

$$\begin{aligned} p_0&= \frac{{{\hat{p}}}_0^m}{tm_{0,x}}, \end{aligned}$$

(29)

$$\begin{aligned} p_\alpha ^q&= \sum _{j=0}^{\alpha }\frac{(-t)^{{\alpha -j}}{{\hat{p}}}_j^q}{z^{\alpha -j+1}}- \frac{{{\hat{p}}}_{0}^m}{m_{0,x}}\sum _{j=0}^{\alpha }\frac{(-t)^{\alpha -j}q_{j,x}}{z^{\alpha -j+1}}, \end{aligned}$$

(30)

$$\begin{aligned} p_\alpha ^m&=\frac{m_{\alpha +1,x}}{m_{0,x}}{{\hat{p}}}_0^m-{{\hat{p}}}_{\alpha +1}^m, \end{aligned}$$

(31)

$$\begin{aligned} p_\beta ^r&= -\sum _{j=0}^{\beta -1}\frac{z^{\beta -j-1}{{\hat{p}}}_j^r}{t^{\beta -j}}+\frac{{{\hat{p}}}_0^m}{m_{0,x}}\left( -\frac{z^{\beta }}{t^{\beta +1}}+\sum _{j=1}^{\beta }\frac{z^{\beta -j}r_{j-1,x}}{t^{\beta -j+1}}\right) , \end{aligned}$$

(32)

where $\alpha \in {\mathbb {N}}_0, \beta \in {\mathbb {N}}$, and ${{\hat{p}}}_0^r=-{{\hat{p}}}_0$.

Moreover, both $({{{\mathcal {R}}}}^{q})^{-1}$, resp. $({{{\mathcal {R}}}}^{m})^{-1}$ leave the set $\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ invariant, i.e., they are recursion operators for $\tau ^W$-symmetries of the rYME (1).

Proof

To prove that ${{{\mathcal {R}}}}^q$, resp. ${{{\mathcal {R}}}}^m$ are bijective, it is enough to verify that $({{{\mathcal {R}}}}^q)^{-1}$, resp $({{{\mathcal {R}}}}^m)^{-1}$ are both the left and right inverses of ${{{\mathcal {R}}}}^q$, resp. ${{{\mathcal {R}}}}^m$ on ${\mathcal {A}}^{{\mathbb {N}}}$. The proofs of the assertions that $({{{\mathcal {R}}}}^{q})^{-1}$, resp. $({{{\mathcal {R}}}}^{m})^{-1}$ leave the set $\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ invariant can be performed similarly as the proofs of Propositions 1 and 2.

Since all results can be readily verified by straightforward but tedious computations we omit the details here. $\square$

Corollary 2

The recursion operators ${\mathcal {R}}^{q}:\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})\rightarrow \textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ and ${\mathcal {R}}^{m}:\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})\rightarrow \textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ for $\tau ^W$-symmetries of the rYME (1) are automorphisms of the real vector space $\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$. The inverses $({\mathcal {R}}^q)^{-1}$, resp. $({\mathcal {R}}^m)^{-1}$ are determined by the formulas (25)–(28), resp. (29)–(32).

Remark 10

Of course, the inverse recursion operators $({\mathcal {R}}^q)^{-1}$ and $({\mathcal {R}}^m)^{-1}$ can also be viewed as infinite-dimensional matrices, see Appendix C. Moreover, it turns out that if we apply $({\mathcal {R}}^q)^{-1}$, resp. $({\mathcal {R}}^m)^{-1}$ on an arbitrary $\tau ^W$-symmetry with a nontrivial local shadow then the shadow of the image, if nonlocal, can depends only on local variables and nonlocal variables from the covering $\tau ^q$, reps. from the covering $\tau ^m$. Our notation of the recursion operators is based just on this observation.

From the point of view of the action of the presented recursion operators on $\tau ^W$-symmetries of the rYME (1), it is also useful to know whether they pair-wise commute. The answer is as follows.

Proposition 4

The ${\mathcal {A}}$-module endomorphisms ${{{\mathcal {R}}}}^q,{{{\mathcal {R}}}}^m,({{{\mathcal {R}}}}^q)^{-1}$ and $({{{\mathcal {R}}}}^m)^{-1}$ pair-wise commute.

Proof

The operators ${{{\mathcal {R}}}}^q$ and $({{{\mathcal {R}}}}^q)^{-1}$ commute as they are mutually inverse. The same is true for the operators ${{{\mathcal {R}}}}^m$ and $({{{\mathcal {R}}}}^m)^{-1}$.

To show that the operators ${{{\mathcal {R}}}}^q$ and ${\mathcal {R}}^{m}$ commute, we apply their commutator on an arbitrary $P\in {{\mathcal {A}}}^{{\mathbb {N}}}$ and use the ${\mathcal {A}}$-linearity of ${\mathcal {R}}^q$. We obtain

$$\begin{aligned} {[}{{{\mathcal {R}}}}^m,{{{\mathcal {R}}}}^q](P)&=[z{\mathcal {I}}+t{{{\mathcal {R}}}}^q,{{{\mathcal {R}}}}^q](P)=z{\mathcal {I}}({{{\mathcal {R}}}}^q(P))+t{{{\mathcal {R}}}}^q({{{\mathcal {R}}}}^q(P))-{{{\mathcal {R}}}}^q((z{\mathcal {I}}+t{{{\mathcal {R}}}}^q)(P))\\&=z{{{\mathcal {R}}}}^q(P)+t{{{\mathcal {R}}}}^q({{{\mathcal {R}}}}^q(P))-{{{\mathcal {R}}}}^q(zP)-{{{\mathcal {R}}}}^q(t{{{\mathcal {R}}}}^q(P))=0. \end{aligned}$$

The remaining commutativity relations between $({\mathcal R}^m)^{-1}$ and ${\mathcal {R}}^q$, $({{\mathcal {R}}}^q)^{-1}$ and ${\mathcal {R}}^m$, and $({{\mathcal {R}}}^q)^{-1}$ and $({\mathcal {R}}^m)^{-1}$ successively follow from the equalities

$$\begin{aligned} {\mathcal {R}}^m\circ [({{\mathcal {R}}}^m)^{-1}, {\mathcal {R}}^q]&={\mathcal {R}}^m \circ \left( ({{\mathcal {R}}}^m)^{-1} \circ {\mathcal {R}}^q- {{{\mathcal {R}}}}^q\circ ({{{\mathcal {R}}}}^m)^{-1}\right) \\&={\mathcal {R}}^q - {\mathcal {R}}^m \circ {\mathcal {R}}^q \circ ({{\mathcal {R}}}^m)^{-1}\\&= {\mathcal {R}}^q - {\mathcal {R}}^q \circ {\mathcal {R}}^m \circ ({\mathcal {R}}^m)^{-1}={\mathcal {R}}^q-{\mathcal {R}}^q =0,\\ {{{\mathcal {R}}}}^q\circ [({{{\mathcal {R}}}}^q)^{-1},{{{\mathcal {R}}}}^m]&={{{\mathcal {R}}}}^q\circ \left( ({{{\mathcal {R}}}}^q)^{-1}\circ {{{\mathcal {R}}}}^m-{{{\mathcal {R}}}}^m\circ ({{{\mathcal {R}}}}^q)^{-1}\right) \\&={{{\mathcal {R}}}}^m-{{{\mathcal {R}}}}^q\circ {{{\mathcal {R}}}}^m\circ ({{{\mathcal {R}}}}^q)^{-1}\\&={{{\mathcal {R}}}}^m-{{{\mathcal {R}}}}^m\circ {{{\mathcal {R}}}}^q\circ ({{{\mathcal {R}}}}^q)^{-1}={{{\mathcal {R}}}}^m-{{{\mathcal {R}}}}^m = 0,\\ {{{\mathcal {R}}}}^q\circ [({{{\mathcal {R}}}}^q)^{-1},({{{\mathcal {R}}}}^m)^{-1}]&={{{\mathcal {R}}}}^q\circ \left( ({{{\mathcal {R}}}}^q)^{-1}\circ ({{{\mathcal {R}}}}^m)^{-1}-({{{\mathcal {R}}}}^m)^{-1}\circ ({{{\mathcal {R}}}}^q)^{-1}\right) \\&=({{{\mathcal {R}}}}^m)^{-1}-{{{\mathcal {R}}}}^q\circ ({{{\mathcal {R}}}}^m)^{-1}\circ ({{{\mathcal {R}}}}^q)^{-1}\\&=({{{\mathcal {R}}}}^m)^{-1}-({{{\mathcal {R}}}}^m)^{-1}\circ {{{\mathcal {R}}}}^q\circ ({{{\mathcal {R}}}}^q)^{-1}\\&=({{{\mathcal {R}}}}^m)^{-1}-({{{\mathcal {R}}}}^m)^{-1} = 0, \end{aligned}$$

and from the fact that ${\mathcal {R}}^q$ and ${\mathcal {R}}^m$ are injective. $\square$

Hereafter, we will call the recursion operators ${\mathcal {R}}^q$, ${\mathcal {R}}^m$ and their inversions $({\mathcal {R}}^q)^{-1}$, $({\mathcal {R}}^m)^{-1}$ as the basic recursion operators; the reason will become clear at the end of this paragraph. The compositions of the basic recursion operators will be denoted as

$$\begin{aligned} {\mathcal {R}}_i^j = \underbrace{({{\mathcal {R}}}^q)^{-1} \circ \ldots \circ ({{\mathcal {R}}}^q)^{-1}}_{i-times} \circ \underbrace{{\mathcal {R}}^m \circ \ldots \circ {\mathcal {R}}^m}_{j-times}, \qquad i,j \in {\mathbb {Z}}, \end{aligned}$$

where for negative values of i, resp. j, we define

$$\begin{aligned}&\underbrace{({{\mathcal {R}}}^q)^{-1} \circ \ldots \circ ({{\mathcal {R}}}^q)^{-1}}_{i-times}=\underbrace{{\mathcal {R}}^q \circ \ldots \circ {\mathcal {R}}^q}_{-i-times}, \\ \mathrm {resp.}\ {}&\underbrace{{\mathcal {R}}^m \circ \ldots \circ {\mathcal {R}}^m}_{j-times}=\underbrace{({{\mathcal {R}}}^m)^{-1} \circ \ldots \circ ({{\mathcal {R}}}^m)^{-1}}_{-j-times}. \end{aligned}$$

In this notation

$$\begin{aligned} {\mathcal {R}}^q \equiv {\mathcal {R}}_{-1}^0, \qquad ({{\mathcal {R}}}^q)^{-1} \equiv {\mathcal {R}}_1^0, \qquad {\mathcal {R}}^m \equiv {\mathcal {R}}_0^1, \qquad ({{\mathcal {R}}}^m)^{-1} \equiv {\mathcal {R}}_0^{-1},\quad {\mathcal {I}} \equiv {\mathcal {R}}_0^0. \end{aligned}$$

It can be easily deduced from Proposition 4 that all the just defined recursion operators ${\mathcal {R}}_i^j,\ i,j\in {\mathbb {Z}}$, pair-wise commute, and, for all $i,j,k,l \in {\mathbb {Z}}$ we have

$$\begin{aligned} {\mathcal {R}}_i^j \circ {\mathcal {R}}_k^l = {\mathcal {R}}_k^l \circ {\mathcal {R}}_i^j = {\mathcal {R}}_{i+k}^{j+l} \mathrm {\ and\ } ({\mathcal {R}}_i^j)^{-1}={\mathcal {R}}_{-i}^{-j}. \end{aligned}$$

Thus, the set ${\{{\mathcal {R}}_i^j\}}_{i,j\in {\mathbb {Z}}}$ carries a commutative group structure with respect to the composition operation, the operators ${\mathcal {R}}_0^1$ and ${\mathcal {R}}_1^0$ forming the minimal set of its generators.

Moreover, it seems that all the recursion operators for $\tau ^W$-symmetries of the rYME (1) arise just as ${\mathbb {R}}$-linear combinations of the operators from $\{ {\mathcal {R}}_i^j \}_{i,j \in {\mathbb {Z}}}$. Thus, we conjecture that the recursion operators ${\mathcal {R}}_{-1}^0, {\mathcal {R}}_1^0, {\mathcal {R}}_0^1$ and ${\mathcal {R}}_0^{-1}$ form the minimal set of generators of the whole associative ${\mathbb {R}}$-algebra of recursion operators for $\tau ^W$-symmetries (with respect to the choice ${\mathcal {A}}=\overline{{\mathcal {L}}\cup {\mathcal {F}}}$) of the rYME (1).

Finally, note that using the defining formula (22) for ${\mathcal {R}}_0^1\equiv {\mathcal {R}}^m$, one can simply obtain an alternative description of the operators ${\mathcal {R}}_i^j$:

$$\begin{aligned} {\mathcal {R}}_i^j = {\mathcal {R}}_0^1 \circ {\mathcal {R}}_i^{j-1} = (t {\mathcal {R}}_{-1}^0 + z {\mathcal {R}}_0^0) \circ {\mathcal {R}}_i^{j-1} = t {\mathcal {R}}_{i-1}^{j-1} + z {\mathcal {R}}_i^{j-1}. \end{aligned}$$

5 Action of the Recursion Operators on $\tau ^W$-Symmetries of the rYME

Once we have the opportunity to grasp the recursion operator as a usual set-theoretical mapping, we could start studying the general properties of its action on the Lie algebra $\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ of $\tau ^W$-symmetries of the rYME (1), as well as we could take a particular symmetry $P \in \textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$ and study its orbit with respect to the action of the group $\{{\mathcal {R}}_{i}^{j} \}_{i,j \in {\mathbb {Z}}}$ on the Lie algebra $\textrm{sym}_{{\mathcal {A}}}^{\tau ^W}({\mathcal {E}})$, which is nothing but a hierarchy of symmetries related to the seed symmetry P. We could then study relations between symmetries within one hierarchy as well as relations between symmetries belonging to distinct hierarchies. However, we intend to address this extensive topic in a separate article.

The principal goal of this section is mainly to point out the computational advantages of using full-fledged recursion operators for the construction of hierarchies of nonlocal symmetries compared to the traditionally used recursion operator for shadows and to demonstrate other benefits and possibilities associated with the existence of a full-fledged operator. Since for the rYME (1) we already have at our disposal both the traditional operators ${\mathcal {R}}^q_{sh}$ and ${\mathcal {R}}^m_{sh}$ for shadows (see Eqs. (21) and (24)) and full-fledged recursion operators, it is natural to use this equation as a basis for such a comparison.

From the computational point of view, it is reasonable to take, if possible, a $\tau ^W$-symmetry with a local shadow as the seed symmetry for each hierarchy. We use two such $\tau ^W$-symmetries as the basis of the first three examples below. Further, let us note that the rYME (1) admits also $\tau ^W$-symmetries involving arbitrary functions $A=A(z,q_0)$, $B=B(t,y)$ and $C=C(z/t,m_0)$ (for brevity we refer to the families of these symmetries as $\textrm{Fam}_A,$ $\textrm{Fam}_B,$ and $\textrm{Fam}_C$). The exact description of these families is a rather complex matter and we also postpone this task to the next paper. Nevertheless, we consider selected symmetries from these families in the last three examples below and we demonstrate how the full-fledged recursion operators may help us with this description.

To achieve the comparison set out above, we first recall the traditional route how one could produce a hierarchy of $\tau ^W$-symmetries for the rYME (1) using the recursion operators ${\mathcal {R}}^q_{sh}$ and ${\mathcal {R}}^m_{sh}$ for shadows.

• (1) We solve the Eq. (NS1) with respect to the unknown (local, if it exists) differential function $p_0$. In other words, we find a (local) shadow $p_0$ of a $\tau ^W$-symmetry of the rYME (1).

• (2) We apply one of the recursion operators for shadows ${\mathcal {R}}^q_{sh}$ and ${\mathcal {R}}^m_{sh}$ to the known shadow $p_0$ to produce a new shadow ${{\hat{p}}}_0$. However, applying any of them to a known shadow $p_0$ and trying to find a new shadow ${{\hat{p}}}_0$, we face to a system of two differential equations, and to find its solution can be in general a very non-trivial task, since the solution may depend on nonlocal variables as well.

• (3) We try to lift each shadow $p_0$ at hand to the full-fledged $\tau ^W$-symmetry, i.e., having fixed a function $p_0$ that satisfies the condition (NS1), we try to solve Eqs. (NS2)–(NS9) with respect to unknown functions $p_{\alpha }^q,\ p_{\alpha }^m$ and $p_{\beta }^r$, $\alpha \ge 0,\ \beta \ge 1.$ Note that the lift must be looked for for each known shadow separately, and moreover, given a shadow, we do not know in advance whether its $\tau ^W$-lift exists at all.

Let us illustrate by an example the disadvantages of this standard approach.

Example 1

By solving Eq. (NS1), one can arrive at the local shadow

$$\begin{aligned} \omega _{0,0}^0=tu_t-xu_x+2u \end{aligned}$$

of a $\tau ^W$-symmetry of the rYME (1). Subsequently, we apply, for example, the recursion operator ${\mathcal {R}}_{sh}^q$ on $\omega _{0,0}^0$ in order to obtain a new shadow. This application can be done in the following two mutually inverse ways.

• (a) The first possibility is to put $p_0=\omega _{0,0}^0$ into Eqs. (21). The new shadow ${{\hat{p}}}_0 = {\mathcal {R}}_{sh}^q(\omega _{0,0}^0)$, which appears in this direction, is determined by the system of differential equations

  $$\begin{aligned} \begin{aligned} {{\hat{p}}}_{0,x}&= (tu_t+2u)u_{xx}-tu_xu_{tx}-xu_{xy}+tu_{ty}-u_x^2+2u_y,\\ {{\hat{p}}}_{0,z}&= xu_zu_{xx}-(tu_z+x)u_{tx}+(tu_t-xu_x+2u)u_{xz}+tu_{tt}-u_xu_z+3u_t. \end{aligned} \end{aligned}$$

  (33)

  To solve such system of equations (recall that we solve it on the covering space $\tilde{{\mathcal {E}}}^W$) is obviously a tough task, since the solution may also depend on nonlocal variables. Indeed, it turns out that the general solution of the system (33) is

  $$\begin{aligned} {{\hat{p}}}_0 = 3r_1+2uu_x+t(u_tu_x+r_{1,t})-xu_y+F(t,y), \end{aligned}$$

  where F(t, y) is an arbitrary function which can be interpreted as an image of the trivial shadow. Thus, here we immediately see another unpleasant feature of the recursion operators for shadows, that they often produce ambiguous results.

• (b) The second possible way is to substitute ${{\hat{p}}}_0=\omega _{0,0}^0$ into Eqs. (21) and search for a new shadow $p_0=({\mathcal {R}}_{sh}^q)^{-1}(\omega _{0,0}^0)$ which appears as a solution of the system of differential equations

  $$\begin{aligned} \begin{aligned} u_xp_{0,x}-p_{0,y}-u_{xx}p_0&= xu_{xx}-tu_{tx}-u_x,\\ u_zp_{0,x}-p_{0,t}-u_{xz}p_0&= xu_{xz}-tu_{tz}-2u_z. \end{aligned} \end{aligned}$$

  (34)

  So, once again we have to deal with the difficult task of finding a general solution, which in this case reads

  $$\begin{aligned} p_0 = t u_z-x + \frac{G(z,q_0)}{q_{0,x}}, \end{aligned}$$

  where $G(z,q_0)$ is an arbitrary function, and the term $G(z,q_0)/q_{0,x}$ arises as an image of the trivial shadow in this direction. Thus, the result is again ambiguous.

Problems of a similar kind are also encountered in the situation when we try to lift the obtained shadows to full-fledged $\tau ^W$-symmetries. For example, we need to do this if we want to identify some hierarchies of commuting $\tau ^W$-symmetries since it is impossible to compute the Jacobi brackets of the $\tau ^W$-symmetries from their shadows only. This step requires solving a system of Eqs. (NS2)–(NS9) for each shadow separately, which in turn involves performing a huge amount of annoying calculations. Moreover, the results of the lifts are again ambiguous since they are given up to the addition of invisible $\tau ^W$-symmetries.

Thus, we can see that to obtain a new shadow, let alone a full-fledged symmetry, from a known shadow using the recursion operator for shadows, we have to perform a highly nontrivial chain of calculations. Moreover, if we want to construct a hierarchy of shadows or full-fledged symmetries, we have to repeat the whole procedure over and over again. However, once we have full-fledged recursion operators at our disposal, we are able to obtain the hierarchy of $\tau ^W$-symmetries of the rYME (1) in a significantly simpler and much more straightforward way: we only once solve Eqs. (NS1)–(NS9) in order to obtain a seed $\tau ^W$-symmetry; other new $\tau ^W$-symmetries (thus, also their shadows) are then generated simply by an iterative application of our basic operators ${\mathcal {R}}_1^0,\ {\mathcal {R}}_0^1,\ {\mathcal {R}}_{-1}^0$ and ${\mathcal {R}}_0^{-1}$.

The wonderful simplicity of the action of the full-fledged recursion operators for symmetries of the rYME (1) is illustrated in the two examples below. In the first one, we present the action of the basic full-fledged recursion operators on the seed full-fledged symmetry $\Omega _0^0$ whose shadow $\omega _{0,0}^0=tu_t-xu_x+2u$ has already been discussed in Ex. 1. In Example 3 below, we illustrate the simplicity of the construction of new full-fledged $\tau ^W$-symmetries of the rYME (1) (and thus also their shadows) via the action of the basic full-fledged recursion operators on another full-fledged $\tau ^W$-symmetry with a local shadow.

Example 2

Upon putting $p_0=\omega _{0,0}^0=tu_t-xu_x+2u$ into (NS2)–(NS9) and solving the resulting equations for unknown functions $p_{\alpha }^q = \omega _{0,\alpha }^{0,q}$, $p_{\alpha }^m = \omega _{0,\alpha }^{0,m}$ and $p_{\beta }^r = \omega _{0,\beta }^{0,r}$, the shadow $\omega _{0,0}^0$ can be lifted to the full-fledged $\tau ^W$-symmetry of the rYME (1)

$$\begin{aligned} \Omega _0^0 =\left[ \omega _{0,0}^0,\omega _{0,0}^{0,q},\omega _{0,0}^{0,m},\omega _{0,1}^{0,r}, \omega _{0,1}^{0,q},\omega _{0,1}^{0,m},\omega _{0,2}^{0,r}, \ldots \right] , \end{aligned}$$

with components

$$\begin{aligned} \omega _{0,\alpha }^{0,q}&= tq_{\alpha ,t}-xq_{\alpha ,x}-\alpha q_\alpha ,\qquad \omega _{0,\alpha }^{0,m} = tm_{\alpha ,t}-xm_{\alpha ,x},\\ \omega _{0,\beta }^{0,r}&=tr_{\beta ,t}-xr_{\beta ,x}+(\beta +2)r_\beta , \qquad \alpha \ge 0, \beta \ge 1. \end{aligned}$$

Let us point out that we omit all possible invisible symmetries of the rYME (1) which could be added to $\Omega _0^0$ in this step, cf. the ending note in Ex. 1.

The simple application of ${\mathcal {R}}_1^0$ provides us with the symmetry $\Omega _1^0= {\mathcal {R}}_1^0(\Omega _0^0)$ with the shadow

$$\begin{aligned} \omega _{1,0}^0 = ({\mathcal {R}}_{1}^0(\Omega _0^0))_1 = \frac{\omega _{0,0}^{0,q}}{q_{0,x}} = \frac{tq_{0,t}-xq_{0,x}}{q_{0,x}} = \frac{tu_zq_{0,x}-xq_{0,x}}{q_{0,x}} = tu_z-x, \end{aligned}$$

and the further components $\omega _{1,\alpha }^{0,q},$ $\omega _{1,\alpha }^{0,m},$ $\omega _{1,\alpha +1}^{0,r}$, $\alpha \ge 0$, in the form:

$$\begin{aligned} \omega _{1,\alpha }^{0,q}&= ({\mathcal {R}}_{1}^0(\Omega _0^0))_{3\alpha +2}\\ {}&=-\frac{q_{\alpha +1,x}}{q_{0,x}}(tq_{0,t}-xq_{0,x}) + tq_{\alpha +1,t}-xq_{\alpha +1,x}-(\alpha +1) q_{\alpha +1}\\ {}&= -t u_z q_{\alpha +1,x}+t(u_z q_{\alpha +1,x}+q_{\alpha ,z})-(\alpha +1)q_{\alpha +1}\\ {}&= t q_{\alpha ,z} -(\alpha +1)q_{\alpha +1},\\ \omega _{1,\alpha }^{0,m}&= ({\mathcal {R}}_1^0(\Omega _0^0))_{3\alpha +3} \\ {}&= \sum \limits _{j=0}^\alpha \frac{t(tm_{j,t}-xm_{j,x})}{(-z)^{\alpha -j+1}} -\frac{t(tq_{0,t}-xq_{0,x})}{q_{0,x}} \cdot \sum \limits _{j=0}^\alpha \frac{m_{j,x}}{(-z)^{\alpha -j+1}}\\ {}&= \sum \limits _{j=0}^\alpha \frac{t^2(m_{j,t}-u_zm_{j,x})}{(-z)^{\alpha -j+1}} = \sum \limits _{j=0}^\alpha \frac{-t(zm_{j,z}+m_{j-1,z})}{(-z)^{\alpha -j+1}}\\ {}&= t \cdot \sum \limits _{j=0}^\alpha \frac{m_{j,z}}{(-z)^{\alpha -j}} - t \cdot \sum \limits _{j=0}^\alpha \frac{m_{j-1,z}}{(-z)^{\alpha -j+1}} = t m_{\alpha ,z},\\ \omega _{1,\alpha +1}^{0,r}&= ({\mathcal {R}}_1^0(\Omega _0^0))_{3\alpha +4} = \frac{r_{\alpha ,x}}{q_{0,x}}(tq_{0,t}-xq_{0,x})-(tr_{\alpha ,t}-xr_{\alpha ,x}+(\alpha +2)r_\alpha )\\ {}&= tu_zr_{\alpha ,x}-tr_{\alpha ,t} - (\alpha +2)r_\alpha = tr_{\alpha +1,z}-(\alpha +2)r_\alpha . \end{aligned}$$

Subsequently, we can very easily obtain another symmetry with a local shadow, namely

$$\begin{aligned} \Omega _1^1&={\mathcal {R}}_0^1(\Omega _1^0)=(t{\mathcal {R}}_{-1}^0+z{\mathcal {R}}_0^0)(\Omega _1^0)=t\, \Omega _0^0+z\,\Omega _1^0, \end{aligned}$$

its shadow being $\omega _{1,0}^1= t\omega _{0,0}^0+z{\omega _{1,0}^0} = t(tu_t-xu_x+2u)+z(tu_z-x)$.

It can be shown that the iterative application of the basic recursion operators ${\mathcal {R}}_{-1}^0,{\mathcal {R}}_{1}^0,{\mathcal {R}}_{0}^1$ and ${\mathcal {R}}_{0}^{-1}$ on $\Omega _0^0$ results in the infinite hierarchy ${\{\Omega _k^l\}}_{k,l \in {\mathbb {Z}}}$ of $\tau ^W$-symmetries of the rYME (1), where $\Omega _k^l={\mathcal {R}}_k^l(\Omega _0^0)$. In Table 1, we present nonlocal shadows of selected $\tau ^W$-symmetries from this hierarchy.

Table 1 Shadows of selected $\tau ^W$-symmetries from the hierarchy $\{\Omega _k^l\}_{k,l \in {\mathbb {Z}}}$

Example 3

One can simply verify that the rYME (1) admits the shadow of symmetry $\psi _{0,0}^0 = u_z$. Upon putting $p_0= \psi _{0,0}^0 =u_z$ into (NS2)–(NS9) and solving the resulting equations for unknown functions $p_{\alpha }^q = \psi _{0,\alpha }^{0,q}$, $p_{\alpha }^m = \psi _{0,\alpha }^{0,m}$ and $p_{\beta }^r = \psi _{0,\beta }^{0,r}$, the shadow $\psi _{0,0}^0$ can be lifted to the full-fledged $\tau ^W$-symmetry of the rYME (1)

$$\begin{aligned} \Psi _0^0 = [u_z, \ldots , q_{\alpha ,z}, m_{\alpha ,z}-(\alpha +1)m_{\alpha +1}, r_{\alpha +1,z}, \ldots ], \quad \alpha \ge 0. \end{aligned}$$

The application of the operator ${\mathcal {R}}_{-1}^0$ to the $\tau ^W$-symmetry $\Psi _0^0$ results in the symmetry

$$\begin{aligned} \Psi _{-1}^0&=\left[ u_t, \ldots , q_{\alpha ,t}, m_{\alpha ,t}+\frac{\alpha }{t} m_\alpha +\frac{(\alpha +1)z}{t}m_{\alpha +1}, r_{\alpha +1,t}, \ldots \right] , \; \alpha \ge 0, \end{aligned}$$

its components being computed as follows:

$$\begin{aligned} \psi _{-1,0}^0&= ({\mathcal {R}}_{-1}^0 (\Psi _0^0))_1 = u_x u_z+r_{1,z} = u_x u_z - u_x u_z + u_t = u_t,\\ \psi _{-1,\alpha }^{0,q}&= ({\mathcal {R}}_{-1}^0 (\Psi _0^0))_{3\alpha +2}=q_{\alpha ,x}u_z+q_{\alpha -1,z} = q_{\alpha ,t},\\ \psi _{-1,\alpha }^{0,m}&= ({\mathcal {R}}_{-1}^0 (\Psi _0^0))_{3\alpha +3}=m_{\alpha ,x}u_z-\frac{1}{t} (m_{\alpha -1,z}-\alpha m_{\alpha })\\&\quad -\frac{z}{t}(m_{\alpha ,z}-(\alpha +1)m_{\alpha +1})\\&= m_{\alpha ,x}u_z - \frac{1}{t} m_{\alpha -1,z} - \frac{z}{t} m_{\alpha ,z} + \frac{\alpha }{t}m_{\alpha } + \frac{(\alpha +1)z}{t}m_{\alpha +1}\\&= m_{\alpha ,t} + \frac{\alpha }{t}m_{\alpha } + \frac{(\alpha +1)z}{t}m_{\alpha +1},\\ \psi _{-1,\alpha +1}^{0,r}&= ({\mathcal {R}}_{-1}^0 (\Psi _0^0))_{3\alpha +4} = r_{\alpha +1,x}u_z-r_{\alpha +2,z} \\&= r_{\alpha +1,x}u_z-r_{\alpha +1,x}u_z + r_{\alpha +1,t} = r_{\alpha +1,t}. \end{aligned}$$

Similarly to Ex. 2, it can be shown again that the iterative application of the basic recursion operators ${\mathcal {R}}_{-1}^0,{\mathcal {R}}_{1}^0,{\mathcal {R}}_{0}^1$ and ${\mathcal {R}}_{0}^{-1}$ on $\Psi _0^0$ results in the infinite hierarchy ${\{\Psi _k^l\}}_{k,l \in {\mathbb {Z}}}$ of $\tau ^W$-symmetries of the rYME (1), where $\Psi _k^l={\mathcal {R}}_k^l(\Psi _0^0)$. For instance, the symmetry $\Psi _0^1$ can be computed very simply using the formula

$$\begin{aligned} \Psi _0^1&= R_0^1(\Psi _0^0)=(t{\mathcal {R}}_{-1}^0+z{\mathcal {R}}_0^0)(\Psi _0^0)=t\Psi _{-1}^0+z\Psi _0^0. \end{aligned}$$

We can immediately see that its shadow is the local function $\psi _{0,0}^1= t\psi _{-1,0}^0+z\psi _{0,0}^0 = tu_t+zu_z$. In Table 2, we present nonlocal shadows of selected $\tau ^W$-symmetries from this hierarchy.

Table 2 Shadows of selected $\tau ^W$-symmetries from the hierarchy ${\{\Psi _k^l\}}_{k,l \in {\mathbb {Z}}}$

The full-fledged recursion operators may prove even more valuable in cases we want to describe at least some of the $\tau ^W$-symmetries involving arbitrary functions, i.e., $\tau ^W$-symmetries from the aforementioned families $\textrm{Fam}_A$, $\textrm{Fam}_B$, resp. $\textrm{Fam}_C$. As it is shown in the following three examples, the great benefit of the full-fledged recursion operators lies in the fact that we are able to straightforwardly describe very complicated $\tau ^W$-symmetries of the rYME (1) that arise as images of much simpler seed symmetries with many zero components. In contrast, if we knew only the shadows and tried to arrive at the same formulas by direct lifting of them, we would have to go through lots of tedious calculations; the success would not be guaranteed at all.

Example 4

The direct computations show that one of the $\tau ^W$-symmetries of the rYME (1) belonging to $\textrm{Fam}_B$ mentioned above is

$$\begin{aligned} \Upsilon _0^0(B) = \left[ \upsilon _{0,0}^0(B), 0,0,\upsilon _{0,1}^{0,r}(B), 0,0, \upsilon _{0,2}^{0,r}(B),0,0, \upsilon _{0,3}^{0,r}(B), \ldots \right] , \end{aligned}$$

where $\upsilon _{0,0}^0(B)=B$, and its other nontrivial components $\upsilon _{0,\beta }^{0,r}(B)$, $\beta \ge 1$ can be quite easily described in a concise form by means of the operator

$$\begin{aligned} {\mathfrak {V}} = z\mathchoice{\frac{\partial }{\partial t}}{\partial /\partial t}{\partial /\partial t}{\partial /\partial t}+x\mathchoice{\frac{\partial }{\partial y}}{\partial /\partial y}{\partial /\partial y}{\partial /\partial y}-2u\mathchoice{\frac{\partial }{\partial x}}{\partial /\partial x}{\partial /\partial x}{\partial /\partial x}+3r_1\mathchoice{\frac{\partial }{\partial u}}{\partial /\partial u}{\partial /\partial u}{\partial /\partial u}-\sum \limits _{i=1}^\infty (i+3)r_{i+1}\mathchoice{\frac{\partial }{\partial r_i}}{\partial /\partial r_i}{\partial /\partial r_i}{\partial /\partial r_i}, \end{aligned}$$

which yields the relations

$$\begin{aligned} \upsilon _{0,\beta }^{0,r}(B) =\frac{(-1)^{\beta -1}}{\beta !}{\mathfrak {V}}^\beta (\upsilon _{0,0}^{0}(B))=\frac{(-1)^{\beta -1}}{\beta !}{\mathfrak {V}}^\beta (B), \quad \beta \ge 1, \end{aligned}$$

where ${\mathfrak {V}}^\beta = \underbrace{{\mathfrak {V}} \circ \ldots \circ {\mathfrak {V}}}_{\beta -times}$.

Applying the basic recursion operator ${\mathcal {R}}_{-1}^0$ on $\Upsilon _0^0(B)$, one can directly obtain the $\tau ^W$-symmetry

$$\begin{aligned} \Upsilon _{-1}^0(B)= & {} \left[ \upsilon _{-1,0}^0(B), \upsilon _{-1,0}^{0,q}(B),\upsilon _{-1,0}^{0,m}(B),\upsilon _{-1,1}^{0,r}(B), \right. \\{} & {} \left. \upsilon _{-1,1}^{0,q}(B),\upsilon _{-1,1}^{0,m}(B), \upsilon _{-1,2}^{0,r}(B), \ldots \right] , \end{aligned}$$

where for each $\alpha \ge 0$, $\beta \ge 1$ we have

$$\begin{aligned} \upsilon _{-1,0}^0(B)&= u_x \upsilon _{0,0}^{0}(B) + \upsilon _{0,1}^{0,r}(B) = u_xB + zB_t + xB_y,\\ \upsilon _{-1,\alpha }^{0,q}(B)&= q_{\alpha ,x} \upsilon _{0,0}^{0}(B) = q_{\alpha ,x} B,\\ \upsilon _{-1,\alpha }^{0,m}(B)&= m_{\alpha ,x} \upsilon _{0,0}^{0}(B) = m_{\alpha ,x} B,\\ \end{aligned}$$

and

$$\begin{aligned} \upsilon _{-1,\beta }^{0,r}(B)&= r_{\beta ,x} \upsilon _{0,0}^{0}(B) - \upsilon _{0,\beta +1}^{0,r}(B) = r_{\beta ,x} B - \frac{(-1)^\beta }{(\beta +1)!} {\mathfrak {V}}^{\beta +1}(B). \end{aligned}$$

Let us remark that both the first column of the matrix ${\mathcal {R}}^q$ (see Sect. 2) and the first column of the matrix ${\mathcal {R}}^m$ (see Appendix C) coincide with the symmetry $\Upsilon _{-1}^0(B)$, the first one upon putting $B(t,y) \equiv 1$, the second one upon putting $B(t,y)=t$.

Example 5

One of the $\tau ^W$-symmetries of the rYME (1) lying in the family $\textrm{Fam}_A$ is the following invisible $\tau ^W$-symmetry:

$$\begin{aligned} \Xi _0^0(A) = \left[ 0, \xi _{0,0}^{0,q}(A),0,0, \xi _{0,1}^{0,q}(A),0,0,\xi _{0,2}^{0,q}(A),0,0, \ldots \right] , \end{aligned}$$

where $\xi _{0,0}^{0,q}(A)=A$, and the components $\xi _{0,\alpha }^{0,q}(A)$, $\alpha \ge 1$, are described by the formula

$$\begin{aligned} \xi _{0,\alpha }^{0,q}(A) = \frac{1}{\alpha }{\mathfrak {X}}(\xi _{0,\alpha -1}^{0,q}(A)) = \frac{1}{\alpha !}{\mathfrak {X}}^\alpha (\xi _{0,0}^{0,q}(A)) = \frac{1}{\alpha !}{\mathfrak {X}}^\alpha (A), \end{aligned}$$

where

$$\begin{aligned} {\mathfrak {X}} = t\mathchoice{\frac{\partial }{\partial z}}{\partial /\partial z}{\partial /\partial z}{\partial /\partial z}+\sum \limits _{i=0}^{\infty } (i+1)q_{i+1}\mathchoice{\frac{\partial }{\partial q_i}}{\partial /\partial q_i}{\partial /\partial q_i}{\partial /\partial q_i}. \end{aligned}$$

The action of the basic recursion operator ${\mathcal {R}}_1^0$ on $\Xi _0^0(A)$ leads to the symmetry

$$\begin{aligned} \Xi _1^0(A) = \left[ \xi _{1,0}^0(A), \xi _{1,0}^{0,q}(A),\xi _{1,0}^{0,m}(A),\xi _{1,1}^{0,r}(A), \xi _{1,1}^{0,q}(A),\xi _{1,1}^{0,m}(A), \xi _{1,2}^{0,r}(A), \ldots \right] , \end{aligned}$$

where for each $\alpha \ge 0$, $\beta \ge 1$ we have

$$\begin{aligned} \xi _{1,0}^{0}(A)&= \frac{\xi _{0,0}^{0,q}(A)}{q_{0,x}}=\frac{A}{q_{0,x}},\\ \xi _{1,\alpha }^{0,q}(A)&= -\frac{q_{\alpha +1,x}}{q_{0,x}}\cdot \xi _{0,0}^{0,q}(A) + \xi _{0,\alpha +1}^{0,q}(A)= -\frac{q_{\alpha +1,x}}{q_{0,x}}\cdot A + \frac{1}{(\alpha +1)!}{\mathfrak {X}}^{\alpha +1}(A),\\ \xi _{1,\alpha }^{0,m}(A)&=-\frac{t \cdot \xi _{0,0}^{0,q}(A)}{q_{0,x}} \cdot \sum \limits _{j=0}^\alpha \frac{m_{j,x}}{(-z)^{\alpha -j+1}}=-\frac{t A}{q_{0,x}} \cdot \sum \limits _{j=0}^\alpha \frac{m_{j,x}}{(-z)^{\alpha -j+1}}\\ \xi _{1,\beta }^{0,q}(A)&=\frac{r_{\beta -1,x}}{q_{0,x}}\cdot \xi _{0,0}^{0,q}(A) =\frac{r_{\beta -1,x}}{q_{0,x}}\cdot A. \end{aligned}$$

Note that upon putting $A(z,q_0) \equiv 1$, we obtain the $\tau ^W$-symmetry $\Xi _{1}^0(1)$ which coincides with the second column of the matrix $({\mathcal {R}}^q)^{-1}$ (see Appendix C).

Example 6

An example of a $\tau ^W$-symmetry of the rYME (1) depending on an arbitrary function $C=C(z/t,m_0)$ is an invisible $\tau ^W$-symmetry

$$\begin{aligned} \Theta _0^0(C) = \left[ 0,0,\theta _{0,0}^{0,m}(C), 0,0, \theta _{0,0}^{0,m}(C),0,0, \theta _{0,0}^{0,m}(C), 0,0,\ldots \right] , \end{aligned}$$

where $\theta _{0,0}^{0,m}(C)=C$, and the components $\theta _{0,\alpha }^{0,m}(C)$, $\alpha \ge 1$ are given by the formula

$$\begin{aligned} \theta _{0,\alpha }^{0,m}(C) = \frac{1}{\alpha }{\mathfrak {T}}(\theta _{0,\alpha -1}^{0,m}(C)) = \frac{1}{\alpha !}{\mathfrak {T}}^\alpha (\theta _{0,0}^{0,m}(C)) = \frac{1}{\alpha !}{\mathfrak {T}}^\alpha (C), \end{aligned}$$

where

$$\begin{aligned} {\mathfrak {T}} = \mathchoice{\frac{\partial }{\partial z}}{\partial /\partial z}{\partial /\partial z}{\partial /\partial z}+\sum \limits _{i=0}^\infty (i+1) m_{i+1}\mathchoice{\frac{\partial }{\partial m_i}}{\partial /\partial m_i}{\partial /\partial m_i}{\partial /\partial m_i}. \end{aligned}$$

Applying the basic recursion operator ${\mathcal {R}}_0^{-1}$ on $\Theta _0^0(C)$, we obtain the $\tau ^W$-symmetry

$$\begin{aligned} \Theta _0^{-1}(C)= & {} \left[ \theta _{0,0}^{-1}(C), \theta _{0,0}^{-1,q}(C),\theta _{0,0}^{-1,m}(C),\theta _{0,1}^{-1,r}(C), \theta _{0,1}^{-1,q}(C),\right. \\{} & {} \left. \theta _{0,1}^{-1,m}(C), \theta _{0,2}^{-1,r}(C), \ldots \right] , \end{aligned}$$

its components for all $\alpha \ge 0$, $\beta \ge 1$ being given by the formulas

$$\begin{aligned} \theta _{0,0}^{-1}(C)&= \frac{\theta _{0,0}^{0,m}(C)}{tm_{0,x}} = \frac{C}{tm_{0,x}},\\ \theta _{0,\alpha }^{-1,q}(C)&= - \frac{\theta _{0,0}^{0,m}(C)}{m_{0,x}}\sum _{j=0}^{\alpha }\frac{(-t)^{\alpha -j}q_{j,x}}{z^{\alpha -j+1}} = - \frac{C}{m_{0,x}}\sum _{j=0}^{\alpha }\frac{(-t)^{\alpha -j}q_{j,x}}{z^{\alpha -j+1}},\\ \theta _{0,\alpha }^{-1,m}(C)&= \frac{m_{\alpha +1,x}}{m_{0,x}} \theta _{0,0}^{0,m}(C) -\theta _{0,\alpha +1}^{0,m}(C) = \frac{m_{\alpha +1,x}}{m_{0,x}} \cdot C -\frac{1}{(\alpha +1)!} {\mathfrak {T}}^{\alpha +1}(C),\\ \end{aligned}$$

and

$$\begin{aligned} \theta _{0,\beta }^{-1,r}(C)&= \frac{\theta _{0,0}^{0,m}(C)}{m_{0,x}}\left( -\frac{z^{\beta }}{t^{\beta +1}}+\sum _{j=1}^{\beta }\frac{z^{\beta -j}r_{j-1,x}}{t^{\beta -j+1}}\right) \\&= \frac{C}{m_{0,x}}\left( -\frac{z^{\beta }}{t^{\beta +1}}+\sum _{j=1}^{\beta }\frac{z^{\beta -j}r_{j-1,x}}{t^{\beta -j+1}}\right) . \end{aligned}$$

It appears that the $\tau ^W$-symmetry $\Theta _0^{-1}(1)$ that lies in the just discussed hierarchy corresponding to the function $C(t/z,m_0)\equiv 1$ coincides with the third column of the matrix $({\mathcal {R}}^m)^{-1}$ (see Appendix C).

6 Conclusions and Future Prospects

In the current literature, recursion operators for nonlocal symmetries of integrable systems are traditionally presented as relations between their shadows only. In principle, knowing such a recursion operator and one shadow, we are able to generate (infinite) hierarchy of new shadows. However, given a differential covering $\tau$ of the equation in question, in case we are interested not only in the shadows but in the full-fledged $\tau$-nonlocal symmetries, we have to investigate separately for each shadow whether it can be lifted to a full-fledged $\tau$-nonlocal symmetry or not. Moreover, another disadvantage of the recursion operators for shadows is that one shadow can be related to more than one (possibly to infinitely many) shadows, since the trivial shadow can be (and usually is) related to a nontrivial one. From this point of view, recursion operators for shadows can be understood to provide us with a ‘highly multivalued mapping’ between full-fledged symmetries. Thus, the difficulties in using such operators to construct an infinite hierarchy of full-fledged $\tau$-nonlocal symmetries become evident.

In this paper, we have shown that the expansion of the Lax pair of a given Lax-integrable equation may provide us with a differential covering which subsequently enables us to find a mapping that puts full-fledged (in this covering) nonlocal symmetries to other ones. This full-fledged recursion operator allows us to generate a (possibly infinite) hierarchy of full-fledged nonlocal symmetries from one known seed full-fledged symmetry directly, without any problem with lifting of shadows. Even though the idea of how to find a full-fledged recursion operator is quite straightforward, the implementation of the latter is not simple at all. The main difficulty of this construction consists in the fact that the covering we have used is infinite-dimensional, and the generating functions of symmetries are in fact infinite sequences of nonlocal functions. Since the recursion operator we are looking for is to act on an infinite sequence resulting in another infinite sequence, our success depends on whether we are able to find general formulas describing the relationships between all the components of the pre-image and its image. We do not claim that this problem can always be satisfactorily overcome; however, we have demonstrated on a particular example of the reduced quasi-classical self-dual Yang–Mills equation that at least sometimes it is definitely possible.

As we have shown in Sects. 2 and 3, the rYME (1) admits two basic pairs of mutually inverse full-fledged recursion operators, such that all of them pair-wise commute. Moreover, it is possible to represent all these recursion operators by infinite-dimensional matrices of differential functions where each matrix has only finitely many nonzero entries in each of its rows, the action of such recursion operators on full-fledged nonlocal symmetries being then simply given as matrix multiplications without any convergence troubles. This enables us to easily generate infinite hierarchies of full-fledged nonlocal symmetries from known ones, see Sect. 4. The use of a full-fledged recursion operator turns out to be particularly advantageous when we let it act on full-fledged symmetries depending on arbitrary functions and generate hierarchies thereof, see Examples 4–6. Nevertheless, a more detailed description of such hierarchies is beyond the scope of this paper and we postpone it, as well as the investigation of the Lie algebra structure on the set of $\tau ^W$-symmetries of the rYME (1) and a more detailed study of the action of the full-fledged recursion operators on it for a forthcoming separate paper.

Based on our preliminary computations, we are pretty sure that the full-fledged recursion operators can be similarly constructed and described also for other linearly degenerate (in the sense of [13]) Lax-integrable equations. As already mentioned in the introductory section of this paper, no serious problems in this direction are expected for equations from Tab. 2 in [9].

Note that recursion operators for full-fledged $\tau$-nonlocal symmetries of the rYME (1) similar to ${\mathcal {R}}^q$ and ${\mathcal {R}}^m$ could be in fact constructed also in the ‘poorer’ covering $\tau =\tau ^q\oplus \tau ^r$. However, the lack of nonlocal variables would make the mapping ${\mathcal {R}}^m$ not surjective, thus, the hierarchies of nonlocal symmetries obtained from one seed symmetry would not be so large. This situation occurs for example in the case we consider the 4D universal hierarchy equation

$$\begin{aligned} u_{yz} = u_{tt}-u_z u_{tx}+u_t u_{xz} \end{aligned}$$

from Tab. 3 in [9], cf. also [7, 15], since it is unclear to us at the moment, given its multi-parameter Lax pair, how one could derive a covering which would be analogous to the $\tau ^m$-covering (15). Let us also remark that it would be interesting to find out whether the construction of a full-fledged recursion operator similar to the one presented here can be done for an equation that is not linearly degenerate, for example for one of the equations that belong to the class of equations of the Monge–Ampère type, see, e.g., [8, 9] and the references therein for further details and examples of such equations.

Finally, recall that one of the frequently investigated properties of recursion operators is their Nijenhuis torsion which is defined, roughly speaking, through the action of a given recursion operator on arbitrary symmetries and their Lie brackets. Hence, verifying the Nijenhuis property of a recursion operator that acts only on shadows is a very tricky task since, in order to compute the Lie brackets, we need to know the full-fledged forms of nonlocal symmetries. However, in the case when the recursion operator for nonlocal symmetries is understood as a linear mapping acting on generating functions of (local) symmetries of a covering equation, the situation is different: we work with the full-fledged forms of these general symmetries and, theoretically, we should be able to express directly their Lie brackets as well as the images of these Lie brackets with respect to the action of the recursion operator in question. Thus, it seems that the full-fledged forms of recursion operators could provide us with the most straightforward way how to compute the Nijenhuis torsion of the latter. We plan to address this issue in more detail in future works.

Appendices

Appendix A: Proof of Proposition 1

We need to check that all components ${{\hat{p}}}_0, {{\hat{p}}}_{\alpha }^q, {{\hat{p}}}_{\alpha }^m, {{\hat{p}}}_{\beta }^r,$ $\alpha \in {\mathbb {N}}_0,$ $\beta \in ~{\mathbb {N}}$, of ${{\hat{P}}}$ satisfy the corresponding conditions (NS1)–(NS9) on the equation $\tilde{{\mathcal {E}}}^W$, provided that $p_0,p_{\alpha }^q,p_{\alpha }^m,p_{\beta }^r$ satisfy them.

• (1) One obtains by straightforward computation:

  $$\begin{aligned} {{\tilde{D}}}_{yz}({{\hat{p}}}_0)&-{{\tilde{D}}}_{tx}({{\hat{p}}}_0) +u_z \tilde{D}_{xx}({{\hat{p}}}_0)-u_x{{\tilde{D}}}_{xz}({{\hat{p}}}_0)+u_{xx}{{\tilde{D}}}_z({{\hat{p}}}_0)-u_{xz}{{\tilde{D}}}_x({{\hat{p}}}_0)\\&=p_0{\tilde{D}}_x(\underbrace{u_{yz}-u_{tx}+u_zu_{xx}-u_xu_{xz}}_{=\, 0 \ \mathrm {due\ to\ Eq.}(1)})\\&\quad +({{\tilde{D}}}_y-u_x{{\tilde{D}}}_x+u_{xx})(\underbrace{\tilde{D}_z(p_1^r)-{{\tilde{D}}}_t(p_0)+u_x{{\tilde{D}}}_z(p_0)+u_z\tilde{D}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS7)})\\&\quad -({{\tilde{D}}}_t-u_z\tilde{D}_x+u_{xz})(\underbrace{{\tilde{D}}_x(p_1^r)-\tilde{D}_y(p_0)+2u_x{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS6)}) =0; \end{aligned}$$

  thus, ${{\hat{p}}}_0$ satisfies the condition (NS1) on the equation $\tilde{{\mathcal {E}}}^W$. This means that ${{\hat{p}}}_0$ is a shadow of nonlocal symmetry od the rYME (1).

• (2) For each $\alpha \in {\mathbb {N}}_0$, we have (recall that $q_j = p_j^q=0$ if $j<0$):

  $$\begin{aligned} {{\tilde{D}}}_t({{\hat{p}}}_{\alpha }^q)&-u_z{{\tilde{D}}}_x({{\hat{p}}}_{\alpha }^q)-q_{\alpha ,x}{{\tilde{D}}}_z({{\hat{p}}}_0)-{{\tilde{D}}}_z({{\hat{p}}}_{\alpha -1}^q)\\&=p_0\tilde{D}_x(\underbrace{q_{\alpha ,t}-u_zq_{\alpha ,x}-q_{\alpha -1,z}}_{=\, 0 \ \mathrm {due\ to\ Eq.}(14)})\\&\quad -q_{\alpha ,x}(\underbrace{{{\tilde{D}}}_z(p_1^r)-\tilde{D}_t(p_0)+u_x{{\tilde{D}}}_z(p_0)+u_z{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS7)})\\&\quad +(\underbrace{{{\tilde{D}}}_t(p_{\alpha -1}^q)-u_z\tilde{D}_x(p_{\alpha -1}^q)-q_{\alpha -1,x}{{\tilde{D}}}_z(p_0)-\tilde{D}_z(p_{\alpha -2}^q)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS2)})=0, \end{aligned}$$

  and also

  $$\begin{aligned} {{\tilde{D}}}_y({{\hat{p}}}_{\alpha }^q)&-u_x{{\tilde{D}}}_x({{\hat{p}}}_{\alpha }^q)-q_{\alpha ,x}{{\tilde{D}}}_x({{\hat{p}}}_0)-{{\tilde{D}}}_x({{\hat{p}}}_{\alpha -1}^q)\\&=p_0\tilde{D}_x(\underbrace{q_{\alpha ,y}-u_xq_{\alpha ,x}-q_{\alpha -1,x}}_{=\, 0\ \mathrm {due\ to\ Eq.}(14)})\\&\quad -q_{\alpha ,x}(\underbrace{{{\tilde{D}}}_x(p_1^r)-\tilde{D}_y(p_0)+2u_x{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS6)})\\&\quad +(\underbrace{{{\tilde{D}}}_y(p_{\alpha -1}^q)-u_x\tilde{D}_x(p_{\alpha -1}^q)-q_{\alpha -1,x}{{\tilde{D}}}_x(p_0)-\tilde{D}_x(p_{\alpha -2}^q)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS3)})=0; \end{aligned}$$

  thus, ${{\hat{p}}}_\alpha ^q$ satisfies the conditions (NS2) and (NS3) on the equation $\tilde{{\mathcal {E}}}^W$.

• (3) For each $\alpha \in {\mathbb {N}}_0$, we have (recall that $m_j = p_j^m=0$ if $j<0$):

  $$\begin{aligned} {{\tilde{D}}}_t({{\hat{p}}}_{\alpha }^m)&-u_z{{\tilde{D}}}_x({{\hat{p}}}_{\alpha }^m)-m_{\alpha ,x}{{\tilde{D}}}_z({{\hat{p}}}_0)+\frac{z}{t}\tilde{D}_z({{\hat{p}}}_{\alpha }^m)+\frac{1}{t}{{\tilde{D}}}_z({{\hat{p}}}_{\alpha -1}^m)\\&=p_0{{\tilde{D}}}_x\underbrace{\left( m_{\alpha ,t}-u_zm_{\alpha ,x}+\frac{z}{t}m_{\alpha ,z}+\frac{1}{t}m_{\alpha -1,z}\right) }_{=\, 0 \ \mathrm {due\ to\ Eq.}(15)}\\&\quad -m_{\alpha ,x} (\underbrace{{{\tilde{D}}}_z(p_1^r)-\tilde{D}_t(p_0)+u_x{{\tilde{D}}}_z(p_0)+u_z{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS7)})\\&\quad -\frac{1}{t} \underbrace{\left( {{\tilde{D}}}_t(p_{\alpha -1}^m)-u_z\tilde{D}_x(p_{\alpha -1}^m)-m_{\alpha -1,x}{{\tilde{D}}}_z(p_0)+\frac{z}{t}\tilde{D}_z(p_{\alpha -1}^m)+\frac{1}{t}{{\tilde{D}}}_z(p_{\alpha -2}^m)\right) }_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS4)}\\&\quad -\frac{z}{t} \underbrace{\left( {{\tilde{D}}}_t(p_{\alpha }^m)-u_z\tilde{D}_x(p_{\alpha }^m)-m_{\alpha ,x}{{\tilde{D}}}_z(p_0)+\frac{z}{t}\tilde{D}_z(p_{\alpha }^m)+\frac{1}{t}{{\tilde{D}}}_z(p_{\alpha -1}^m)\right) }_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS4)} = 0, \end{aligned}$$

  and also

  $$\begin{aligned} {{\tilde{D}}}_y({{\hat{p}}}_{\alpha }^m)&-u_x{{\tilde{D}}}_x({{\hat{p}}}_{\alpha }^m)-m_{\alpha ,x}{{\tilde{D}}}_x({{\hat{p}}}_0)+\frac{z}{t}\tilde{D}_x({{\hat{p}}}_{\alpha }^m)+\frac{1}{t}{{\tilde{D}}}_x({{\hat{p}}}_{\alpha -1}^m)\\&=p_0{{\tilde{D}}}_x\underbrace{\left( m_{\alpha ,y}-u_xm_{\alpha ,x}+\frac{z}{t}m_{\alpha ,x}+\frac{1}{t}m_{\alpha -1,x}\right) }_{=\, 0 \ \mathrm {due\ to\ Eq.}(15)}\\&\quad -m_{\alpha ,x} (\underbrace{{{\tilde{D}}}_x(p_1^r)-\tilde{D}_y(p_0)+2u_x{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS6)})\\&\quad -\frac{1}{t} \underbrace{\left( {{\tilde{D}}}_y(p_{\alpha -1}^m)-u_x\tilde{D}_x(p_{\alpha -1}^m)-m_{\alpha -1,x}{{\tilde{D}}}_x(p_0)+\frac{z}{t}\tilde{D}_x(p_{\alpha -1}^m)+\frac{1}{t}{{\tilde{D}}}_x(p_{\alpha -2}^m)\right) }_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS5)}\\&\quad -\frac{z}{t} \underbrace{\left( {{\tilde{D}}}_y(p_{\alpha }^m)-u_x\tilde{D}_x(p_{\alpha }^m)-m_{\alpha ,x}{{\tilde{D}}}_x(p_0)+\frac{z}{t}\tilde{D}_x(p_{\alpha }^m)+\frac{1}{t}{{\tilde{D}}}_x(p_{\alpha -1}^m)\right) }_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS5)} = 0; \end{aligned}$$

  thus, ${{\hat{p}}}_\alpha ^m$ satisfies the conditions (NS4) and (NS5) on the equation $\tilde{{\mathcal {E}}}^W$.

• (4) For the quantity ${{\hat{p}}}_1^r$, we obtain:

  $$\begin{aligned} {{\tilde{D}}}_x({{\hat{p}}}_1^r)&-{{\tilde{D}}}_y({{\hat{p}}}_0)+2u_x{{\tilde{D}}}_x({{\hat{p}}}_0)\\&=p_0{{\tilde{D}}}_x(\underbrace{r_{1,x}+u_x^2-u_y}_{=\, 0 \ \mathrm {due\ to\ Eq.}(16)})+u_x (\underbrace{{{\tilde{D}}}_x(p_1^r)-\tilde{D}_y(p_0)+2u_x{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS6)})\\&\quad -(\underbrace{{{\tilde{D}}}_x(p_2^r)-u_x{{\tilde{D}}}_x(p_1^r)-r_{1,x}\tilde{D}_x(p_0)+{{\tilde{D}}}_y(p_1^r)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS8)}) =0, \end{aligned}$$

  and also

  $$\begin{aligned} {{\tilde{D}}}_z({{\hat{p}}}_1^r)&-{{\tilde{D}}}_t({{\hat{p}}}_0)+u_x{{\tilde{D}}}_z({{\hat{p}}}_0)+u_z{{\tilde{D}}}_x({{\hat{p}}}_0)\\&=p_0{{\tilde{D}}}_x(\underbrace{r_{1,z}+u_xu_z-u_t}_{=\, 0 \ \mathrm {due\ to\ Eq.}(16)})\\&\quad +u_x(\underbrace{{{\tilde{D}}}_z(p_r^1)-{{\tilde{D}}}_t(p_0)+u_x\tilde{D}_z(p_0)+u_z{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS7)})\\&\quad -(\underbrace{{{\tilde{D}}}_z(p_2^r)-u_z{{\tilde{D}}}_x(p_1^r)-r_{1,x}\tilde{D}_z(p_0)+{{\tilde{D}}}_t(p_1^r)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS9)}) = 0; \end{aligned}$$

  thus, ${{\hat{p}}}_1^r$ satisfies the conditions (NS6) and (NS7) on the equation $\tilde{{\mathcal {E}}}^W$.

• (5) For each $\beta \ge 2$, we have:

  $$\begin{aligned} {{\tilde{D}}}_x({{\hat{p}}}_{\beta }^r)&-u_x{{\tilde{D}}}_x({{\hat{p}}}_{\beta -1}^r)-r_{\beta -1,x}{{\tilde{D}}}_x({{\hat{p}}}_0)+{{\tilde{D}}}_y({{\hat{p}}}_{\beta -1}^r)\\&=p_0\tilde{D}_x(\underbrace{r_{\beta ,x}-u_xr_{\beta -1,x}+r_{\beta -1,y}}_{=\, 0 \ \mathrm {due\ to\ Eq.}(16)}) \\&\quad - r_{\beta -1,x}(\underbrace{{{\tilde{D}}}_x(p_1^r)-\tilde{D}_y(p_0)+2u_x{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS6)})\\&\quad -(\underbrace{{{\tilde{D}}}_x(p_{\beta +1}^r)-u_x{{\tilde{D}}}_x( p_{\beta }^r)-r_{\beta ,x}{{\tilde{D}}}_x(p_0)+\tilde{D}_y(p_{\beta }^r)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS8)})=0, \end{aligned}$$

  and also

  $$\begin{aligned} {{\tilde{D}}}_z({{\hat{p}}}_{\beta }^r)&-u_z{{\tilde{D}}}_x({{\hat{p}}}_{\beta -1}^r)-r_{\beta -1,x}{{\tilde{D}}}_z({{\hat{p}}}_0)+{{\tilde{D}}}_t({{\hat{p}}}_{\beta -1}^r)\\&=p_0{{\tilde{D}}}_x (\underbrace{r_{\beta ,z}-u_zr_{\beta -1,x}+r_{\beta -1,t}}_{=\, 0 \ \mathrm {due\ to\ Eq.}(16)})\\&\quad -r_{\beta -1,x}(\underbrace{{{\tilde{D}}}_z(p_1^r)-\tilde{D}_t(p_0)+u_x{{\tilde{D}}}_z(p_0)+u_z{{\tilde{D}}}_x(p_0)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS7)})\\&\quad -(\underbrace{{{\tilde{D}}}_z(p_{\beta +1}^r)-u_z{{\tilde{D}}}_x( p_{\beta }^r)-r_{\beta ,x}{{\tilde{D}}}_z(p_0)+\tilde{D}_t(p_{\beta }^r)}_{=\, 0 \ \mathrm {due\ to\ Eq.}(NS9)}) =0; \end{aligned}$$

  thus, ${{\hat{p}}}_\beta ^r$ satisfies the conditions (NS8) and (NS9) on the equation $\tilde{{\mathcal {E}}}^W$.

$\square$

Appendix B: Proof of Proposition 2

Let $P \in {\textrm{sym}}_{{\mathcal {L}}}^{\tau ^W}({\mathcal {E}})$ be an arbitrary $\tau ^W$-symmetry of the rYME (1), and denote ${{\hat{P}}} = {\mathcal {R}}^q(P)$, ${{\bar{P}}} = {\mathcal {R}}^m(P)$. Then ${{\hat{P}}} \in {\textrm{sym}}_{{\mathcal {L}}}^{\tau ^W}({\mathcal {E}})$ and, according to (22), it holds

$$\begin{aligned} {{\bar{P}}}={\mathcal {R}}^m (P)=(z{\mathcal {I}}+t{\mathcal {R}}^q)(P)=t {{\hat{P}}} + zP. \end{aligned}$$

We need to show that ${{\bar{P}}} \in {\textrm{sym}}_{{\mathcal {L}}}^{\tau ^W}({\mathcal {E}})$.

For brevity, let us denote by $L_i$, $i=1,\ldots ,9$, the left-hand sides of the respective conditions (NS1)–(NS9). Then $L_i(P)=L_i({{\hat{P}}})=0$ on $\tilde{{\mathcal {E}}}^W$ ($L_i(P)$ means the substitution of the corresponding components of P to $L_i$, etc.), and one can obtain by straightforward computations for each $\alpha \ge 0,\beta \ge 2$ the following:

$$\begin{aligned} L_1({{\bar{P}}})&=tL_1({{\hat{P}}})+zL_1(P)-{\tilde{D}}_x({{\hat{p}}}_0)+{\tilde{D}}_y(p_0)-u_x{\tilde{D}}_x(p_0)+u_{xx} p_0\\&=-{{\tilde{D}}}_x(p_1^r)+{\tilde{D}}_y(p_0)-2u_x{\tilde{D}}_x(p_0) = -L_6(P)=0,\\ L_2({{\bar{P}}})&=tL_2({{\hat{P}}}) + zL_2(P) + \underbrace{{{\hat{p}}}_{\alpha }^q-q_{\alpha ,x}p_0-p_{\alpha -1}^q}_{=0 \ \mathrm {due\ to\ Eq.} (18)} = 0,\\ L_3({{\bar{P}}})&=tL_3({{\hat{P}}}) + zL_3(P) = 0\\ L_4({{\bar{P}}})&=tL_4({{\hat{P}}}) + zL_4(P) +\underbrace{{{\hat{p}}}_{\alpha }^m-m_{\alpha ,x}p_0+\frac{z}{t}p_{\alpha }^m+\frac{1}{t} p_{\alpha -1}^m}_{=0 \ \mathrm {due\ to\ Eq.} (19)} = 0\\ L_5({{\bar{P}}})&=tL_5({{\hat{P}}}) + zL_5(P) = 0\\ L_6({{\bar{P}}})&=tL_6({{\hat{P}}}) + zL_6(P) = 0\\ L_7({{\bar{P}}})&=tL_7({{\hat{P}}}) + zL_7(P) - \underbrace{{{\hat{p}}}_0+p_1^r+u_xp_0}_{=0 \ \mathrm {due\ to\ Eq.} (17)}=0\\ L_8({{\bar{P}}})&=tL_8({{\hat{P}}}) + zL_8(P) = 0\\ L_9({{\bar{P}}})&=tL_9({{\hat{P}}}) + zL_9(P) + \underbrace{{{\hat{p}}}_{\beta -1}^r+p_{\beta }^r-r_{{\beta }-1,x}p_0}_{=0 \ \mathrm {due\ to\ Eq.} (20)}= 0. \end{aligned}$$

where $p_{-1}^q=p_{-1}^m=0$. Thus, the components of ${{\bar{P}}}$ satisfy all conditions (NS1)–(NS9) on $\tilde{{\mathcal {E}}}^W$ which implies ${{\bar{P}}} \in {\textrm{sym}}_{{\mathcal {L}}}^{\tau ^W}({\mathcal {E}})$. $\square$

Appendix C: Left-upper corners of the matrix representations of the recursion operators $({{\mathcal {R}}}^q)^{-1}$, ${\mathcal {R}}^m$, and $({{\mathcal {R}}}^m)^{-1}$

$$\begin{aligned}{} & {} ({{\mathcal {R}}}^q)^{-1}\\{} & {} \quad = \left( \begin{array}{lllllllllllll} 0 &{} 1/q_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} -q_{1,x}/q_{0,x} &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} -\frac{t}{q_{0,x}} \cdot \sum \limits _{j=0}^0 \frac{m_{j,x}}{(-z)^{1-j}} &{} -\frac{t}{z} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 1 &{} r_{0,x}/q_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} -q_{2,x}/q_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} -\frac{t}{q_{0,x}} \cdot \sum \limits _{j=0}^1 \frac{m_{j,x}}{(-z)^{2-j}} &{} \frac{t}{z^2} &{} 0 &{} 0 &{} -\frac{t}{z} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} r_{1,x}/q_{0,x} &{} 0 &{} -1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} -q_{3,x}/q_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} \ldots \\ 0 &{} -\frac{t}{q_{0,x}} \cdot \sum \limits _{j=0}^2 \frac{m_{j,x}}{(-z)^{3-j}} &{} -\frac{t}{z^3} &{} 0 &{} 0 &{} \frac{t}{z^2} &{} 0 &{} 0 &{} -\frac{t}{z} &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} r_{2,x}/q_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} -1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \ddots \end{array}\right) \end{aligned}$$

$$\begin{aligned} {{\mathcal {R}}^m}= \left( \begin{array}{lllllllllllllll} z+t u_x &{} 0 &{} 0 &{} t &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tq_{0,x} &{} z &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tm_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tr_{1,x} &{} 0 &{} 0 &{} z &{} 0 &{} 0 &{} -t &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tq_{1,x} &{} t &{} 0 &{} 0 &{} z &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tm_{1,x} &{} 0 &{} -1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tr_{2,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} z &{} 0 &{} 0 &{} -t &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tq_{2,x} &{} 0 &{} 0 &{} 0 &{} t &{} 0 &{} 0 &{} z &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tm_{2,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} -1&{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ tr_{3,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} z &{} 0 &{} 0 &{} -t &{} 0 &{} \ldots \\ \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \ddots \end{array}\right) \end{aligned}$$

$$\begin{aligned}{} & {} ({{\mathcal {R}}}^m)^{-1}\\{} & {} \quad = \left( \begin{array}{llllllllll} 0 &{} 0 &{} 1/(tm_{0,x}) &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} \frac{1}{z} &{} - \frac{1}{m_{0,x}}\sum \limits _{j=0}^{0}\frac{(-t)^{-j}q_{j,x}}{z^{1-j}} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} 0 &{} m_{1,x}/m_{0,x} &{} 0 &{} 0 &{} -1 &{} 0 &{} 0 &{} 0 &{} \ldots \\ \frac{1}{t} &{} 0 &{} \frac{1}{m_{0,x}}\left( -\frac{z}{t^2}+\sum \limits _{j=1}^{1}\frac{z^{1-j}r_{j-1,x}}{t^{2-j}}\right) &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} -\frac{t}{z^2} &{} - \frac{1}{m_{0,x}}\sum \limits _{j=0}^{1}\frac{(-t)^{1-j}q_{j,x}}{z^{2-j}} &{} 0 &{} \frac{1}{z} &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} 0 &{} m_{2,x}/m_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} -1 &{} \ldots \\ \frac{z}{t^2} &{} 0 &{} \frac{1}{m_{0,x}}\left( -\frac{z^2}{t^3}+\sum \limits _{j=1}^{2}\frac{z^{2-j}r_{j-1,x}}{t^{3-j}}\right) &{} -\frac{1}{t} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ 0 &{} \frac{t^2}{z^3} &{}- \frac{1}{m_{0,x}}\sum \limits _{j=0}^{2}\frac{(-t)^{2-j}q_{j,x}}{z^{3-j}} &{} 0 &{} -\frac{t}{z^2} &{} 0 &{} 0 &{} \frac{1}{z} &{} 0 &{} \ldots \\ 0 &{} 0 &{} m_{3,x}/m_{0,x} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \ldots \\ \frac{z^2}{t^3} &{} 0 &{} \frac{1}{m_{0,x}}\left( -\frac{z^3}{t^4}+\sum \limits _{j=1}^{3}\frac{z^{3-j}r_{j-1,x}}{t^{4-j}}\right) &{} -\frac{z}{t^2} &{} 0 &{} 0 &{} -\frac{1}{t} &{} 0 &{} 0 &{} \ldots \\ \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \vdots &{} \ddots \end{array}\right) \end{aligned}$$