text:
# Askey-Wilson type functions, with bound states

## 1. Introduction

Orthogonal polynomials which are eigenfunctions of a differential operator have a long history. When the differential operator is of order two, Bochner (1929)  proved that this property characterizes the so-called classical orthogonal polynomials, linked with the names of Hermite, Laguerre and Jacobi. The general problem was raised by H.L. Krall  in 1938. He proved that the differential operator has to be of even order and in  he obtained a complete classification for the case of an operator of order four.

In , it was observed that Krall’s result can be reformulated in a striking manner. The (semi-infinite) tridiagonal matrices which encode the three-term recursion relation satisfied by Krall’s polynomials, are obtained by adding a bound state (with an arbitrary weight) at the boundary of the continuous spectrum of the tridiagonal matrices associated with some instances of the Laguerre and the Jacobi polynomials. This way of looking at Krall’s result leads to a very efficient construction of his polynomials, by means of the so-called (discrete) Darboux transformation, one of the basic tools in the theory of solitons. Roughly, the method consists in factorizing the tridiagonal matrix (appropriately shifted) as a product of a lower and an upper matrix, producing then a new matrix by permuting the factors. One of the nice outcomes of this approach is that the Darboux process can be iterated, leading to orthogonal sequences of polynomials satisfying differential equations, with a moment functional given by a weight distribution involving not only Dirac’s delta function but also any of its derivatives. We refer the reader to  and  for the precise formulation of the results in the context of the Laguerre and the Jacobi polynomials respectively. For a sample of other works related to Krall’s problem or the role of the Darboux transformation in the theory of orthogonal polynomials, the reader can consult , , , , , ,  and .

This paper deals with a discrete version of Krall’s problem, replacing the differential operator by a $`q`$-difference operator. The polynomials that one meets here are of more recent vintage. The celebrated Askey-Wilson polynomials  form the most general family of orthogonal polynomials which are also eigenfunctions of a second order $`q`$-difference operator, see , . In fact, we address an extension of the problem, by allowing for a *doubly infinite* three-term recursion relation, instead of a semi-infinite one:

###### The Askey-Wilson bispectral problem.

To determine all doubly infinite tridiagonal matrices $`L`$, for which some family of eigenfunctions $`\mathrm{\Psi }(n,z)`$ satisfying

$$A_n\mathrm{\Psi }(n+1,z)+B_n\mathrm{\Psi }(n,z)+C_n\mathrm{\Psi }(n1,z)=(z+z^1)\mathrm{\Psi }(n,z),$$
(1.1)
is also a family of eigenfunctions of a $`q`$-difference operator $`B`$ in the spectral variable $`z`$, i.e.

$$\underset{\text{finitely many }i\text{’s}}{}D_i(z)\mathrm{\Psi }(n,q^iz)=\mathrm{\Lambda }(n)\mathrm{\Psi }(n,z).$$
(1.2)
The motivation for allowing for doubly infinite tridiagonal matrices stems from the pioneering work of Duistermaat and Grünbaum , where a purely differential version of the problem (with the tridiagonal operator $`L`$ replaced by the Schrödinger operator) was posed and completely solved. The problem above was raised in  by one of us in collaboration with F.A. Grünbaum, and solved there in the case of a $`q`$-difference operator of order two. The result is that, by allowing for doubly infinite tridiagonal matrices, the only solutions are provided by what we proposed to call the *associated Askey-Wilson matrices*. The entries of these matrices are obtained by making an *arbitrary* shift in the coefficients of the recursion relations satisfied by the Askey-Wilson polynomials, extending them over *all* integers. Precisely,

$$A_n=\stackrel{~}{A}_{n+\epsilon },B_n=\stackrel{~}{B}_{n+\epsilon },C_n=\stackrel{~}{C}_{n+\epsilon },n,$$
(1.3)
where $`\stackrel{~}{A}_n,\stackrel{~}{B}_n,\stackrel{~}{C}_n`$, denote the coefficients of the three-term recursion relation satisfied by the Askey-Wilson polynomials. The coefficients (1.3) define the so-called *associated Askey-Wilson polynomials* when the conditions $`\mathrm{\Psi }(1,z)=0`$ and $`\mathrm{\Psi }(0,z)=1`$ are imposed in (1.1). The associated Askey-Wilson polynomials were extensively studied by Ismail and Rahman in the remarkable paper . However, as soon as $`\epsilon 0`$, the functions $`\mathrm{\Psi }(\gamma ,z),\gamma n+\epsilon `$, that solve (1.1) and (1.2), with a second-order $`q`$-difference operator, are *not* given by the associated polynomials, but rather by a *two-dimensional* space of contiguous basic hypergeometric *functions*, that we called the *Askey-Wilson functions*. In , these functions were constructed recursively in terms of an arbitrary solution of a second-order $`q`$-difference equation (called there the Gauss-Askey-Wilson equation), but explicit formulas were not obtained.

Our first result, Theorem 2.1, is to give explicit formulas for the Askey-Wilson functions in terms of basic hypergeometric series. For this, we use two linearly independent eigenfunctions of the Askey-Wilson second order difference operator $`LL_{a,b,c,d}(\gamma ,E)`$ in terms of very-well-poised $`{}_{8}{}^{}\varphi _{7}^{}`$ functions, that were found by Ismail and Rahman , . The notation $`L_{a,b,c,d}(\gamma ,E)`$ reminds that the operator $`L`$ depends on four arbitrary parameters (that are usually denoted by $`a,b,c,d`$), $`\gamma `$ stands for $`n+\epsilon `$ and $`E`$ denotes the customary shift operator. We show that by multiplying the Ismail-Rahman functions by an appropriate $`z`$ dependent (but $`\gamma `$ independent) factor, these functions satisfy a beautiful symmetry property, which we shall refer to as duality. It essentially means that the geometric and the spectral parameters $`\gamma `$ and $`z`$ of the eigenfunctions are interchangeable, if we introduce an appropriate involution on the parameters $`a,b,c,d`$. A special case of this duality was considered by Koelink and Stokman , in their study of the Askey-Wilson function transform. It is important to note that these authors focus on a *very specific* eigenfunction of the Askey-Wilson difference operator, which they call the Askey-Wilson function, because it is a meromorphic continuation of the Askey-Wilson polynomial in its degree. The same function was studied previously by Suslov , , who established Fourier-Bessel type orthogonality relations for it.

Sections 3 to 7 form the core of the paper. The main result can be found in Section 6, Theorem 6.2 and Corollary 6.6. We prove that starting with *some* instances of the associated Askey-Wilson matrices, and adding bound states (with an arbitrary weight) at a collection of *specific* points *outside of* the continuous spectrum of these operators, we can construct solutions of the Askey-Wilson bispectral problem, as stated at the beginning of this introduction. We naturally call these solutions *Askey-Wilson type functions*. When $`\epsilon =0`$ in (1.3), we obtain in this way orthogonal polynomials which are eigenfunctions of $`q`$-difference operators of arbitrary orders, thus providing $`q`$-analogues of Krall’s orthogonal polynomials. As shown by formula (6.18) in Theorem 6.2 (see also Corollary 6.6), when $`q1`$, all the added bound states accumulate at the boundary $`\pm 2`$ of the continuous spectrum $`[2,2]`$ of the Askey-Wilson second order difference operator. It is interesting to note that most works on Krall’s original problem deal with measures obtained by adding delta functions (with arbitrary masses) to the measure of orthogonality of some instances of the classical orthogonal polynomials, at the boundary of the interval of orthogonality, discarding the possibility of adding also derivatives of the delta function. To the best of our knowledge, situations involving derivatives of the delta function were first contemplated in . In this case we don’t have an orthogonality *measure*. Our result shows that, in the $`q`$-difference case, the bound states *split up* and we do get an orthogonality measure for the resulting polynomials. We shall now describe our strategy to establish Theorem 6.2, the proof of which is prepared by the material exposed in Sections 3 to 5 of the paper.

The basic technique to establish Theorem 6.2 is the so-called method of bispectral Darboux transformations which was developed by Bakalov, Horozov and Yakimov , , and by Kasman and Rothstein , in relation with a program aiming at describing all bispectral commutative rings of differential operators. In , the method was adapted to attack more systematically Krall’s original problem. Though it was successful to produce extensions of the so-called Krall-Laguerre polynomials, it led quickly to serious computational difficulties, when applied starting with the Jacobi polynomials. In the context of the Jacobi polynomials, the difficulties were overcome by Grünbaum and Yakimov , by introducing a new idea that we shall explain below. Section 3 explains the method of bispectral Darboux transformations. Roughly, the method consists in finding the most general ”rational” factorization of some constant coefficient polynomial in $`L_{a,b,c,d}`$

$$=𝔥(L_{a,b,c,d})=\underset{i=1}{\overset{m}{}}(L_{a,b,c,d}x_i\mathrm{Id})=𝒫𝒬,$$
(1.4)
for appropriate choices of the parameters $`a,b,c,d`$ and of the bound states $`x_i`$ to be added. The factorization involves $`m`$ free parameters defining a new tridiagonal operator $`\widehat{L}`$, by exchanging the two factors

$$=𝔥(L_{a,b,c,d})=𝒫𝒬\widehat{}=𝒬𝒫=𝔥(\widehat{L}).$$

By a ”rational” factorization, it is meant that each of the factors $`𝒫`$ and $`𝒬`$ can be written as a polynomial in the Askey-Wilson difference operator $`L_{a,b,c,d}`$ and the corresponding diagonal operator $`\mathrm{\Lambda }_\gamma `$ in (1.2), after factoring out some rational ”function” in $`\mathrm{\Lambda }_\gamma `$. The precise formulation is to be found in Theorem 3.1. The crucial Theorem 3.3 adapts the Grünbaum-Yakimov technology alluded to above, within the context of our paper. This theorem characterizes the bispectral Darboux transformations described in Theorem 3.1 as those for which each factor $`𝒫`$ and $`𝒬`$ is a difference operator with rational coefficients in $`q^\gamma `$ which is moreover invariant under the involution

$$I(q^\gamma )=\frac{q^{\gamma +1}}{abcd}\text{and}I(E)=E^1.$$
(1.5)
Section 4 expresses some contiguous relations between the Askey-Wilson functions, in the language of the discrete Darboux transformation. This is used later in Section 6 to form in some cases the appropriate polynomial $``$ in (1.4) that needs to be factorized, see Remark 6.4. Section 5 deals with the explicit computation of the kernel of the operator $`𝒬`$, from which the Darboux factorization $`=𝒫𝒬`$ is performed in Section 6. It is quite remarkable that checking the invariance of the operators $`𝒫`$ and $`𝒬`$ under the involution (1.5), involves the use of classical formulas in the theory of basic hypergeometric series, such as the Sears’ and Watson’s transformation formulas. Section 7 illustrates Theorem 6.2 on the simplest possible example, adding one bound state off the continuous spectrum. When $`\epsilon `$ in (1.3) is put to be zero, this example leads to orthogonal polynomials which are eigenfunctions of a $`q`$-difference operator of order 4, providing a $`q`$-deformation of the original Krall-Jacobi polynomials.

Section 8 of the paper requires some background in algebraic geometry. Though it played a decisive role in our research, it is not needed for the understanding of the earlier sections of the paper. Its aim is to provide some further solutions of the Askey-Wilson bispectral problem, as well as to take a first step in putting the problem in the more general context of the duality property satisfied by the Askey-Wilson functions, that we established in Section 2. This requires to allow in the statement of the problem itself for *arbitrary* difference as well as $`q`$-difference operators, leaving thus the context of orthogonal polynomials.

In , we started a systematic study of bispectral commutative rings of difference operators, for which the dual ring is a ring of differential operators, a problem which can be thought of as an extension of Krall’s original question. The techniques we used were adapted from the beautiful paper of G. Wilson , who classified all bispectral rank 1 commutative rings of differential operators. By definition, the rank of a commutative ring of differential or difference operators is the greatest common divisor of the orders of all the operators in the ring. Here again the basic philosophy is ”duality”: rank 1 commutative rings of differential or difference operators enjoy some bispectral property, when the common eigenfunction of the operators in these rings possesses a symmetry which allows for an exchange of the geometric and the spectral variable. Ultimately, this property reflects itself in the geometry of the moduli space of algebraic curves which describes the spectrum (in the sense of algebraic geometry) of these rings. In the context of  and , the spectrum of the corresponding rings must be an affine irreducible *rational* curve with only *cusp-like singularities*.

Applying the techniques of , in Proposition 8.8, we exhibit a class of *rational* curves with *double* points which are the spectrum of rank 1 bispectral commutative rings of difference operators, with dual ring a ring of $`q`$-difference operators. If we further impose to these rings to contain a tridiagonal operator, we are led to the special family of rational curves with double points given by

$$v^2=(u^24)\underset{iJ}{}\left(u(q^{k_i/2}+q^{k_i/2})\right)^2\underset{iK}{}\left(u+(q^{k_i/2}+q^{k_i/2})\right)^2,$$
(1.6)
with $`k_i,1ig`$, an arbitrary collection of positive integers, which have been partitioned in two arbitrary sets $`J`$ and $`K`$. This is the content of Theorem 8.1, which is the central result of Section 8. The corresponding tridiagonal operators depend on $`g`$ arbitrary parameters. They can be obtained by iteration of the Darboux transformation, starting from the operator

$$E+E^1,$$
(1.7)
and adding bound states at the points $`q^{k_i/2}+q^{k_i/2}`$ for $`iJ`$, and $`(q^{k_i/2}+q^{k_i/2})`$ for $`iK`$. It is well known that, when appropriately parametrized, the tridiagonal operators obtained in this way provide solitonic solutions to the doubly infinite Toda lattice hierarchy, see . Thus, we call these solutions *Askey-Wilson type solitons*. They are the analogues of the rational solutions of the Korteweg-de Vries hierarchy, which were found by Duistermaat and Grünbaum  to form one of the two families of solutions to the purely differential version of the bispectral problem. When $`q1`$, the $`q`$-difference bispectral operator (appropriately scaled) in (1.2) becomes a differential operator. All the added bound states accumulate at the boundary $`\pm 2`$ of the continuous spectrum $`[2,2]`$ of the operator (1.7). The corresponding curves (1.6) acquire then only cusp-like singularities. We refer the reader to  for a complete discussion of this limiting case.

In conclusion, the solutions of the Askey-Wilson bispectral problem described in Theorem 8.1, as well as those described in Theorem 6.2 and Corollary 6.6 parallel the so-called rank 1 and rank 2 solutions found by Duistermaat and Grünbaum , for the purely differential version of the problem. Whether we have found all solutions of the Askey-Wilson bispectral problem can only be pure speculation at this point. We feel that finding an appropriate ”moduli space” of basic hypergeometric functions on which our result would follow from duality, looks to be an extremely interesting and challenging problem. Let us only mention that, even within the context of the Askey-Wilson polynomials, the deeper understanding of duality stems from affine Hecke algebraic considerations (see ), which were pioneered by Cherednik ,  in his proof of some conjectures about the Macdonald polynomials.

## 2. Askey-Wilson functions

Throughout the paper we use the standard notations for basic hypergeometric series, following the book  by Gasper and Rahman. In particular, we write

$$(a_1,a_2,\mathrm{},a_r;q)_k=\underset{i=1}{\overset{r}{}}(a_i;q)_k,$$

with

$$(a;q)_k=\frac{(a;q)_{\mathrm{}}}{(aq^k;q)_{\mathrm{}}}\text{and}(a;q)_{\mathrm{}}=\underset{i=0}{\overset{\mathrm{}}{}}(1aq^i),$$

for products of $`q`$-shifted factorials, where $`0<q<1`$. The series expansion

$${}_{r}{}^{}\varphi _{s}^{}[\begin{array}{c}a_1,a_2,\mathrm{},a_r\\ b_1,b_2,\mathrm{},b_s\end{array};q,z]=\underset{k=0}{\overset{\mathrm{}}{}}\frac{(a_1,a_2,\mathrm{},a_r;q)_k}{(q,b_1,\mathrm{},b_s;q)_k}\left[(1)^kq^{\left(\genfrac{}{}{0pt}{}{k}{2}\right)}\right]^{1+sr}z^k$$

defines the $`{}_{r}{}^{}\varphi _{s}^{}`$ basic hypergeometric series.

We denote by $`E`$ and $`E^1`$, respectively, the customary forward and backward shift operators, acting on a function $`h_\gamma =h(\gamma )`$ by

$$Eh_\gamma =h_{\gamma +1},E^1h_\gamma =h_{\gamma 1}.$$

The Askey-Wilson difference operator $`L_{a,b,c,d}(\gamma ,E)`$ is the second-order difference operator

$$L_{a,b,c,d}(\gamma ,E)=A_\gamma E+B_\gamma \mathrm{Id}+C_\gamma E^1,$$
(2.1)
where $`\mathrm{Id}`$ is the identity operator and the coefficients $`A_\gamma =A_{\gamma ;a,b,c,d}`$, $`B_\gamma =B_{a,b,c,d}`$, $`C_\gamma =C_{\gamma ;a,b,c,d}`$ are given by

$`A_\gamma `$ $`={\displaystyle \frac{(1abq^\gamma )(1acq^\gamma )(1adq^\gamma )(1abcdq^{\gamma 1})}{a(1abcdq^{2\gamma 1})(1abcdq^{2\gamma })}},`$ (2.2a)
$`C_\gamma `$ $`={\displaystyle \frac{a(1q^\gamma )(1bcq^{\gamma 1})(1bdq^{\gamma 1})(1cdq^{\gamma 1})}{(1abcdq^{2\gamma 2})(1abcdq^{2\gamma 1})}},`$ (2.2b)
$`B_\gamma `$ $`=a+a^1(A_\gamma +C_\gamma ).`$ (2.2c)
Two linearly independent solutions of the functional equation

$$L_{a,b,c,d}(\gamma ,E)h_\gamma (z)=(z+1/z)h_\gamma (z)$$
(2.3)
in terms of very-well-poised $`{}_{8}{}^{}\varphi _{7}^{}`$ series, were constructed by Ismail and Rahman in , . They are given by

$`r_\gamma (a,b,c,d;z)`$ $`={\displaystyle \frac{(abq^\gamma ,acq^\gamma ,adq^\gamma ,bcdq^\gamma /z;q)_{\mathrm{}}}{(bcq^\gamma ,bdq^\gamma ,cdq^\gamma ,azq^\gamma ;q)_{\mathrm{}}}}\left({\displaystyle \frac{a}{z}}\right)^\gamma \times `$
$`{}_{8}{}^{}W_{7}^{}(bcd/qz;b/z,c/z,d/z,abcdq^{\gamma 1},q^\gamma ;q,qz/a),`$ (2.4)
and
$`s_\gamma (a,b,c,d;z)`$ $`={\displaystyle \frac{(abcdq^{2\gamma },bzq^{\gamma +1},czq^{\gamma +1},dzq^{\gamma +1},bcdzq^\gamma ;q)_{\mathrm{}}}{(bcq^\gamma ,bdq^\gamma ,cdq^\gamma ,q^{\gamma +1},bcdzq^{2\gamma +1};q)_{\mathrm{}}}}(az)^\gamma \times `$
$`{}_{8}{}^{}W_{7}^{}(bcdzq^{2\gamma };bcq^\gamma ,bdq^\gamma ,cdq^\gamma ,q^{\gamma +1},zq/a;q,az),`$ (2.5)
where $`{}_{8}{}^{}W_{7}^{}(a;a_1,a_2,a_3,a_4,a_5;q,z)`$ stands for the very-well-poised hypergeometric series

$${}_{8}{}^{}\varphi _{7}^{}[\begin{array}{c}a,qa^{1/2},qa^{1/2},a_1,a_2,a_3,a_4,a_5\\ a^{1/2},a^{1/2},qa/a_1,qa/a_2,qa/a_3,qa/a_4,qa/a_5\end{array};q,z].$$

Let us denote by $`D_z`$ and $`D_z^1`$, respectively, the forward and backward $`q`$-shift operators, acting on a function $`h(z)`$ by

$$D_zh(z)=h(qz)\text{and}D_z^1h(z)=h(z/q).$$

Finally, let $`B_{a,b,c,d}(z,D_z)`$ be the Askey-Wilson second-order $`q`$-difference operator

$$B_{a,b,c,d}(z,D_z)=A(z)D_z\left[A(z)+A(1/z)\right]\mathrm{Id}+A(1/z)D_z^1,$$
(2.6)
where

$$A(z)=A_{a,b,c,d}(z)=\frac{(1az)(1bz)(1cz)(1dz)}{(1z^2)(1qz^2)}.$$
(2.7)
With these notations we have

###### Theorem 2.1.

The functions

$`R_\gamma (a,b,c,d;z)`$ $`={\displaystyle \frac{(az;q)_{\mathrm{}}}{(bcd/z;q)_{\mathrm{}}}}r_\gamma (a,b,c,d;z)`$ (2.8a)
and
$`S_\gamma (a,b,c,d;z)`$ $`={\displaystyle \frac{z^{1+\mathrm{log}(abcd)/\mathrm{log}(q)}(az;q)_{\mathrm{}}}{(zq/a,zq/b,zq/c,zq/d;q)_{\mathrm{}}}}s_\gamma (a,b,c,d;z)`$ (2.8b)
solve the bispectral problem

$`L_{a,b,c,d}(\gamma ,E)\mathrm{\Psi }(\gamma ,z)`$ $`=(z+1/z)\mathrm{\Psi }(\gamma ,z),`$ (2.9)
$`B_{a,b,c,d}(z,D_z)\mathrm{\Psi }(\gamma ,z)`$ $`=\mathrm{\Lambda }_\gamma \mathrm{\Psi }(\gamma ,z),`$ (2.10)
where

$$\mathrm{\Lambda }_\gamma =q^\gamma (1q^\gamma )(1abcdq^{\gamma 1}).$$
(2.11)
Note that $`R_\gamma (a,b,c,d;z)`$ and $`S_\gamma (a,b,c,d;z)`$ differ from $`r_\gamma (a,b,c,d;z)`$ and $`s_\gamma (a,b,c,d;z)`$, respectively, by a factor independent of $`\gamma `$, hence (2.9) is automatically satisfied. The proof of (2.10) is based on an interesting bispectral involution. As we shall see, when appropriately normalized, $`r_\gamma (a,b,c,d;z)`$ and $`s_\gamma (a,b,c,d;z)`$ satisfy an important duality relation. All this is the content of the next lemma.

Let us introduce “dual” parameters and variables via the formulas

$$\stackrel{~}{a}=\sqrt{\frac{q^3}{abcd}},\stackrel{~}{b}=\frac{q^2}{\stackrel{~}{a}ab},\stackrel{~}{c}=\frac{q^2}{\stackrel{~}{a}ac},\stackrel{~}{d}=\frac{q^2}{\stackrel{~}{a}ad},$$
(2.12a)
$$\stackrel{~}{z}=\frac{q^{\gamma +1}}{\stackrel{~}{a}},\stackrel{~}{\gamma }=\frac{\mathrm{log}(az)}{\mathrm{log}(q)}1,$$
(2.12b)
and define

$`\overline{r}(a,b,c,d;z)`$ $`=bcdz^{2\mathrm{log}(abcd)/\mathrm{log}(q)}{\displaystyle \frac{(bc,bd,cd,qz/a,qz/b,qz/c,qz/d;q)_{\mathrm{}}}{(bz,cz,dz,bcd/z;q)_{\mathrm{}}}},`$ (2.13)
$`\overline{s}(a,b,c,d;z)`$ $`={\displaystyle \frac{bcdz}{(bz,cz,dz;q)_{\mathrm{}}}}.`$ (2.14)
Denote by $`\overline{R}`$ and $`\overline{S}`$ the functions

$`\overline{R}_\gamma (a,b,c,d;z)`$ $`=\overline{r}(a,b,c,d;z)r_\gamma (a,b,c,d;z),`$ (2.15)
$`\overline{S}_\gamma (a,b,c,d;z)`$ $`=\overline{s}(a,b,c,d;z)s_\gamma (a,b,c,d;z).`$ (2.16)
###### Lemma 2.2.

The functions $`\overline{R}`$ and $`\overline{S}`$ defined above satisfy the duality relations

$`\overline{R}_\gamma (a,b,c,d;z)`$ $`=\overline{R}_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z}),`$ (2.17)
and
$`\overline{S}_\gamma (a,b,c,d;z)`$ $`=\overline{S}_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z}).`$ (2.18)
###### Proof of Lemma 2.2.

Both (2.17) and (2.18) can be proved by using transformation formula \[10, (III.24), p. 243\]. Indeed, replacing $`a`$, $`b`$, $`c`$, $`d`$, $`e`$, $`f`$ by $`bcd/qz`$, $`abcdq^{\gamma 1}`$, $`c/z`$, $`d/z`$, $`b/z`$, $`q^\gamma `$ we have

$`{}_{8}{}^{}W_{7}^{}(bcd/qz;b/z,c/z,d/z,abcdq^{\gamma 1},q^\gamma ;q,qz/a)`$
$`={\displaystyle \frac{(bcd/z,abcdq^{\gamma 1},qb/a,qc/a,qd/a,zq^{\gamma +1}/a;q)_{\mathrm{}}}{(bc,bd,cd,bcdq^\gamma /z,q^{\gamma +2}/a^2,qz/a;q)_{\mathrm{}}}}\times `$ (2.19)
$`{}_{8}{}^{}W_{7}^{}(q^{\gamma +1}/a^2;q^{\gamma +1}/ab,q^{\gamma +1}/ac,q^{\gamma +1}/ad,qz/a,q/az;q,abcdq^{\gamma 1}).`$
The $`{}_{8}{}^{}W_{7}^{}`$ series in the right-hand side of (2.19) can be rewritten in terms of the dual parameters (see (2.12)) as

$${}_{8}{}^{}W_{7}^{}(\stackrel{~}{b}\stackrel{~}{c}\stackrel{~}{d}/q\stackrel{~}{z};\stackrel{~}{b}/\stackrel{~}{z},\stackrel{~}{c}/\stackrel{~}{z},\stackrel{~}{d}/\stackrel{~}{z},\stackrel{~}{a}\stackrel{~}{b}\stackrel{~}{c}\stackrel{~}{d}q^{\stackrel{~}{\gamma }1},q^{\stackrel{~}{\gamma }};q,q\stackrel{~}{z}/\stackrel{~}{a}),$$

which is exactly the $`{}_{8}{}^{}W_{7}^{}`$ factor in $`r_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z})`$. Now taking into account (2.4), (2.12), (2.13) and (2.19), one can easily check that

$$\frac{r_\gamma (a,b,c,d;z)}{r_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z})}=\frac{\overline{r}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z})}{\overline{r}(a,b,c,d,z)},$$

which combined with (2.15) gives (2.17). The proof of (2.18) is similar. First we apply \[10, (III.24), p. 243\] with $`a`$, $`b`$, $`c`$, $`d`$, $`e`$, $`f`$ replaced by $`bcdzq^{2\gamma }`$, $`q^{\gamma +1}`$, $`bcq^\gamma `$, $`bdq^\gamma `$, $`cdq^\gamma `$, $`zq/a`$ to obtain

$`{}_{8}{}^{}W_{7}^{}(bcdzq^{2\gamma };bcq^\gamma ,bdq^\gamma ,cdq^\gamma ,q^{\gamma +1},zq/a;q,az)`$
$`={\displaystyle \frac{(bcdzq^{2\gamma +1},q^{\gamma +1},abczq^\gamma ,abdzq^\gamma ,acdzq^\gamma ,qz^2;q)_{\mathrm{}}}{(bzq^{\gamma +1},czq^{\gamma +1},dzq^{\gamma +1},abcdq^{2\gamma },abcdz^2q^\gamma ,az;q)_{\mathrm{}}}}\times `$ (2.20)
$`{}_{8}{}^{}W_{7}^{}(abcdz^2q^{\gamma 1};az,bz,cz,dz,abcdq^{\gamma 1};q,q^{\gamma +1}).`$
Again the hypergeometric series in the right-hand side can be written in the dual parameters as

$${}_{8}{}^{}W_{7}^{}(\stackrel{~}{b}\stackrel{~}{c}\stackrel{~}{d}\stackrel{~}{z}q^{2\stackrel{~}{\gamma }};\stackrel{~}{b}\stackrel{~}{c}q^{\stackrel{~}{\gamma }},\stackrel{~}{b}\stackrel{~}{d}q^{\stackrel{~}{\gamma }},\stackrel{~}{c}\stackrel{~}{d}q^{\stackrel{~}{\gamma }},q^{\stackrel{~}{\gamma }+1},\stackrel{~}{z}q/\stackrel{~}{a};q,\stackrel{~}{a}\stackrel{~}{z}),$$

and the proof follows by combining (2.5), (2.20), (2.12), (2.14) and (2.16). ∎

###### Proof of Theorem 2.1.

Let $`F=F_\gamma (a,b,c,d;z)`$ be $`\overline{R}_\gamma (a,b,c,d)`$ or $`\overline{S}_\gamma (a,b,c,d)`$. From Lemma 2.2, it follows that $`F`$ satisfies also a difference equation in $`\stackrel{~}{\gamma }`$, i.e. we have

$$L_{\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d}}(\stackrel{~}{\gamma },E_{\stackrel{~}{\gamma }})F_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z})=(\stackrel{~}{z}+1/\stackrel{~}{z})F_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z}).$$
(2.21)
On the other hand, from (2.12b), it is clear that

$$F_{\stackrel{~}{\gamma }\pm 1}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d};\stackrel{~}{z})=D_z^{\pm 1}F_\gamma (a,b,c,d;z).$$
(2.22)
Combining (2.1), (2.2c), (2.21) and (2.22) we see that $`F_\gamma (a,b,c,d;z)`$ satisfies the following $`q`$-difference equation in $`z`$

$$\begin{array}{cc}& \stackrel{~}{a}(\stackrel{~}{z}+1/\stackrel{~}{z}\stackrel{~}{a}1/\stackrel{~}{a})F_\gamma (a,b,c,d;z)\hfill \\ & =\overline{A}(z)\left(F_\gamma (a,b,c,d;zq)F_\gamma (a,b,c,d;z)\right)\hfill \\ & \overline{C}(z)\left(F_\gamma (a,b,c,d;z/q)F_\gamma (a,b,c,d;z)\right),\hfill \end{array}$$
(2.23)
where $`\overline{A}(z)=\stackrel{~}{a}A_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d})`$ and $`\overline{C}(z)=\stackrel{~}{a}C_{\stackrel{~}{\gamma }}(\stackrel{~}{a},\stackrel{~}{b},\stackrel{~}{c},\stackrel{~}{d})`$. Using (2.2a), (2.2b), and (2.12), one can express the coefficients of the $`q`$-difference operator and the eigenvalue in formula (2.23) in terms of the parameters $`a`$, $`b`$, $`c`$, $`d`$ and of the variables $`z`$ and $`\gamma `$ as

$`\overline{A}(z)`$ $`={\displaystyle \frac{(1qz/a)(1qz/b)(1qz/c)(1qz/d)}{(1z^2)(1qz^2)}},`$ (2.24a)
$`\overline{C}(z)`$ $`={\displaystyle \frac{q^3}{abcd}}{\displaystyle \frac{(1az/q)(1bz/q)(1cz/q)(1dz/q)}{(1z^2)(1z^2/q)}},`$ (2.24b)
and

$$\stackrel{~}{a}(\stackrel{~}{z}+1/\stackrel{~}{z}\stackrel{~}{a}1/\stackrel{~}{a})=q^{\gamma +1}+\frac{q^{2\gamma }}{abcd}\frac{q^3}{abcd}1.$$
(2.25)
If we define $`\mathrm{\Phi }(z)`$ to be the function

$$\mathrm{\Phi }(z):=z^{2+\mathrm{log}(abcd)/\mathrm{log}(q)}\frac{(az,bz,cz,dz;q)_{\mathrm{}}}{(zq/a,zq/b,zq/c,zq/d;q)_{\mathrm{}}},$$
(2.26)
then from (2.23), (2.24) and (2.25) it follows that $`G(z):=\mathrm{\Phi }(z)F_\gamma (a,b,c,d;z)`$ satisfies the $`q`$-difference equation

$$\left(q^{\gamma +1}+\frac{q^{2\gamma }}{abcd}\frac{q^3}{abcd}1\right)G(z)=\widehat{A}(z)G(zq)\widehat{B}(z)G(z)+\widehat{C}(z)G(z/q),$$
(2.27)
where

$`\widehat{A}(z)`$ $`={\displaystyle \frac{\mathrm{\Phi }(z)}{\mathrm{\Phi }(zq)}}\overline{A}(z)={\displaystyle \frac{q^2}{abcd}}A(z),`$ (2.28a)
$`\widehat{C}(z)`$ $`={\displaystyle \frac{\mathrm{\Phi }(z)}{\mathrm{\Phi }(z/q)}}\overline{C}(z)={\displaystyle \frac{q^2}{abcd}}A(1/z),`$ (2.28b)
$`\widehat{B}(z)`$ $`=\overline{A}(z)+\overline{C}(z),`$ (2.28c)
with $`A(z)`$ defined in (2.7). A direct computation now shows that

$$\widehat{B}(z)=\frac{q^2}{abcd}\left(A(z)+A(1/z)+q1\frac{abcd}{q}+\frac{abcd}{q^2}\right),$$
(2.29)
which combined with (2.27) and (2.28) proves that $`G(z)`$ satisfies the Askey-Wilson $`q`$-difference equation (2.10). To finish the proof, it is enough to notice that

$`\mathrm{\Phi }(z)\overline{R}_\gamma (a,b,c,d)`$ $`=bcd(bc,bd,cd;q)_{\mathrm{}}R_\gamma (a,b,c,d;z),`$ (2.30)
and
$`\mathrm{\Phi }(z)\overline{S}_\gamma (a,b,c,d)`$ $`=bcdS_\gamma (a,b,c,d;z),`$ (2.31)
i.e. $`\mathrm{\Phi }(z)\overline{R}`$ and $`\mathrm{\Phi }(z)\overline{S}`$ differ from $`R`$ and $`S`$, respectively, by unessential constant factors. ∎

## 3. Bispectral Darboux transformations

We denote by

$$=L_{a,b,c,d},\mathrm{\Lambda }_\gamma ,$$
(3.1)
the algebra of difference operators generated by $`L_{a,b,c,d}`$ and $`\mathrm{\Lambda }_\gamma `$, defined by (2.1) and (2.11), respectively. Similarly

$$^{}=z+z^1,B_{a,b,c,d},$$
(3.2)
will denote the algebra of $`q`$-difference operators generated by the operator of multiplication by $`z+z^1`$ and the operator $`B_{a,b,c,d}`$ defined in (2.6). Formulas (2.9) and (2.10) serve to define an anti-isomorphism

$$𝔟:^{}$$
(3.3)
between these two algebras, i.e. it is given on the generators by

$$𝔟(L_{a,b,c,d})=z+z^1\text{and}𝔟(\mathrm{\Lambda }_\gamma )=B_{a,b,c,d}.$$
(3.4)
Note that the bispectral property (2.9), (2.10) is equivalent to the identity

$$X\mathrm{\Psi }(\gamma ,z)=𝔟(X)\mathrm{\Psi }(\gamma ,z),X.$$
(3.5)
We shall also need the commutative subalgebras (the algebras of “functions”) of $``$ and $`^{}`$ defined by

$$𝒦=\mathrm{\Lambda }_\gamma \text{and}𝒦^{}=z+z^1.$$
(3.6)
The next theorem summarizes the technology of bispectral Darboux transformations, initiated in ,  and . It was adapted and applied to the case of difference operators in . For the convenience of the reader, we include the short proof of this result.

###### Theorem 3.1.

Let $``$ be a constant coefficient polynomial in $`L_{a,b,c,d}`$, which factorizes rationally as

$$=𝒫𝒬,$$
(3.7)
in such a way that

$$𝒫=U\mathrm{\Gamma }^1,𝒬=\mathrm{\Theta }^1V,$$
(3.8)
with $`U,V`$, and $`\mathrm{\Theta },\mathrm{\Gamma }𝒦`$. Then the Darboux transform of $``$ given by

$$\widehat{}=𝒬𝒫$$
(3.9)
is a bispectral operator. More precisely, defining $`f=𝔟()𝒦^{}`$ and $`\widehat{\mathrm{\Psi }}=𝒬\mathrm{\Psi }`$, with $`\mathrm{\Psi }`$ satisfying (2.9) and (2.10), we have

$`\widehat{}\widehat{\mathrm{\Psi }}`$ $`=f\widehat{\mathrm{\Psi }},`$ (3.10)
$`\widehat{B}\widehat{\mathrm{\Psi }}`$ $`=\mathrm{\Theta }\mathrm{\Gamma }\widehat{\mathrm{\Psi }},`$ (3.11)
with

$$\widehat{B}=𝔟(V)𝔟(U)f^1.$$
(3.12)
###### Proof.

From the definitions above and equation (3.5), we obtain

$$\widehat{}\widehat{\mathrm{\Psi }}=𝒬𝒫𝒬\mathrm{\Psi }=𝒬\mathrm{\Psi }=𝒬𝔟()\mathrm{\Psi }=f𝒬\mathrm{\Psi }=f\widehat{\mathrm{\Psi }},$$

which establishes (3.10). Using (3.5) and (3.8), we can write $`\widehat{\mathrm{\Psi }}`$ as

$$\widehat{\mathrm{\Psi }}=𝒬\mathrm{\Psi }=\mathrm{\Theta }^1V\mathrm{\Psi }=\mathrm{\Theta }^1𝔟(V)\mathrm{\Psi }.$$
(3.13)
From (3.7) and (3.8), we deduce that

$$\mathrm{\Theta }\mathrm{\Gamma }=V^1U.$$

Applying the anti-isomorphism $`𝔟`$ to this equation, we obtain

$$𝔟(\mathrm{\Theta }\mathrm{\Gamma })=𝔟(U)f^1𝔟(V).$$
(3.14)
Finally, using (3.5), (3.12), (3.13) and (3.14) we get

$`\mathrm{\Theta }\mathrm{\Gamma }\widehat{\mathrm{\Psi }}`$ $`=\mathrm{\Theta }\mathrm{\Gamma }\mathrm{\Theta }^1𝔟(V)\mathrm{\Psi }=\mathrm{\Theta }^1𝔟(V)\mathrm{\Theta }\mathrm{\Gamma }\mathrm{\Psi }=\mathrm{\Theta }^1𝔟(V)𝔟(\mathrm{\Theta }\mathrm{\Gamma })\mathrm{\Psi }`$
$`=\mathrm{\Theta }^1𝔟(V)𝔟(U)f^1𝔟(V)\mathrm{\Psi }=𝔟(V)𝔟(U)f^1\mathrm{\Theta }^1𝔟(V)\mathrm{\Psi }=\widehat{B}\widehat{\mathrm{\Psi }},`$
which gives (3.11) and completes the proof. ∎

Despite the apparent simplicity of Theorem 3.1, it is a priori very complicated to recognize that an operator admits a rational Darboux factorization as defined by (3.7), (3.8). The rest of the section is concerned with making Theorem 3.1 more effective. Some of the ideas developed below were first introduced in , in the context of the Jacobi polynomials.

Let us denote by $`\{q^\gamma ,E\}`$ the algebra of difference operators of the form

$$T=\underset{j=m_1}{\overset{m_2}{}}h_j(q^\gamma )E^j,$$

with coefficients $`h_j(q^\gamma )`$ rational functions in $`q^\gamma `$. The ordered pair $`[m_1,m_2]`$ is called the support of $`T`$. Consider the involution<sup>1</sup><sup>1</sup>1By involution we mean an automorphism $`I`$, such that $`II=\mathrm{Id}`$. $`I`$ on $`\{q^\gamma ,E\}`$, defined by

$$I(q^\gamma )=\frac{q^{\gamma +1}}{abcd}\text{and}I(E)=E^1.$$
(3.15)
###### Lemma 3.2.

A Laurent polynomial in $`q^\gamma `$

$$p(q^\gamma )=\underset{k=m}{\overset{n}{}}c_kq^{k\gamma },m,n,$$
(3.16)
is a polynomial in $`\mathrm{\Lambda }_\gamma `$ (i.e. it belongs to $`𝒦`$ as defined in (3.6)) if and only if it is $`I`$-invariant.

###### Proof.

A straightforward computation from (2.11) gives

$$I(\mathrm{\Lambda }_\gamma )=\mathrm{\Lambda }_\gamma ,$$

hence every polynomial in $`\mathrm{\Lambda }_\gamma `$ is $`I`$-invariant. Conversely, if $`p(q^\gamma )`$ in (3.16) is $`I`$-invariant, one checks that $`m=n`$ and

$$p(q^\gamma )=\underset{k=0}{\overset{n}{}}c_k\left(q^{k\gamma }+\frac{q^{k(1\gamma )}}{(abcd)^k}\right).$$
(3.17)
Observe from (2.11) that

$$q^\gamma +\frac{q^{1\gamma }}{abcd}=\frac{q}{abcd}\left(\mathrm{\Lambda }_\gamma +1+\frac{abcd}{q}\right).$$

Since

$$\begin{array}{c}q^{k\gamma }+\frac{q^{k(1\gamma )}}{(abcd)^k}=\left(q^\gamma +\frac{q^{1\gamma }}{abcd}\right)\left(q^{(k1)\gamma }+\frac{q^{(k1)(1\gamma )}}{(abcd)^{k1}}\right)\hfill \\ \hfill \frac{q}{abcd}\left(q^{(k2)\gamma }+\frac{q^{(k2)(1\gamma )}}{(abcd)^{k2}}\right),\end{array}$$

by induction, it follows that $`p(q^\gamma )`$ in (3.17) is a polynomial in $`\mathrm{\Lambda }_\gamma `$, i.e. it belongs to $`𝒦`$. ∎

An easy computation, using formulas (2.2a) and (2.2b), gives

$$I(A_{\gamma ;a,b,c,d})=C_{\gamma ;a,b,c,d},$$
(3.18)
which shows that the operator $`L_{a,b,c,d}`$ in (2.1) is $`I`$-invariant. Let us denote by $`\mathrm{\Delta }^I`$ the subalgebra of $`\{q^\gamma ,E\}`$, consisting of $`I`$-invariant operators, i.e.

$$\mathrm{\Delta }^I=\left\{T\{q^\gamma ,E\}\right|I(T)=T\}.$$

Since $`L_{a,b,c,d}\mathrm{\Delta }^I`$ and $`\mathrm{\Lambda }_\gamma \mathrm{\Delta }^I`$, it follows that

$$𝒦\mathrm{\Delta }^I.$$
(3.19)
Let $`\mathrm{\Delta }_m^I`$ be the subset of $`\mathrm{\Delta }^I`$, consisting of $`I`$-invariant operators with support $`[m,m]`$. The next theorem characterizes the bispectral Darboux transformations as those for which the two factors $`𝒫`$ and $`𝒬`$ in (3.7) are $`I`$-invariant operators, with rational coefficients.

###### Theorem 3.3.

The following conditions on an operator $`T\{q^\gamma ,E\}`$ are equivalent:

* The operator $`T`$ is $`I`$-invariant, i.e. $`T\mathrm{\Delta }^I`$;
* $`T`$ has the form $`\mathrm{\Theta }^1V`$, for some operator $`V`$ and some function $`\mathrm{\Theta }𝒦`$;
* $`T`$ has the form $`U\mathrm{\Gamma }^1`$, for some operator $`U`$ and some function $`\mathrm{\Gamma }𝒦`$.

###### Proof.

From (3.19) it is obvious that (ii) implies (i). Below, we shall prove that (i) implies (ii), i.e. if $`T\mathrm{\Delta }^I`$, then $`T=\mathrm{\Theta }^1V`$, for some $`\mathrm{\Theta }𝒦`$ and $`V`$. The proof of the fact that (i) is equivalent to (iii) is similar.

It is clear that $`T`$ must have support $`[m,m]`$, for some $`m`$. If $`m=0`$, we can always assume that $`T`$ is the quotient of two $`I`$-invariant Laurent polynomials in $`q^\gamma `$. Hence, by Lemma 3.2, $`T`$ is a rational function in $`\mathrm{\Lambda }_\gamma `$, so we can take $`V`$ and $`\mathrm{\Theta }`$ in (ii) to be the numerator and the denominator, respectively. If $`m1`$, we shall show that

$$T=\mathrm{\Theta }_m^1V_m\text{mod}\mathrm{\Delta }_{m1}^I,$$
(3.20)
with $`\mathrm{\Theta }_m𝒦`$, $`V_m`$, from which the proof follows by induction.

Any operator in $`\mathrm{\Delta }_m^I`$ may be written

$$T=(\text{Id}+I)\frac{𝔠(q^\gamma )}{𝔡(q^\gamma )}E^m\text{mod}\mathrm{\Delta }_{m1}^I,$$
(3.21)
with $`𝔠(q^\gamma )`$ and $`𝔡(q^\gamma )`$ polynomials in $`q^\gamma `$. In particular, since $`L_{a,b,c,d}`$ is $`I`$-invariant, we have

$$L_{a,b,c,d}^m=(\text{Id}+I)\frac{𝔣(q^\gamma )}{𝔤(q^\gamma )}E^m\text{mod}\mathrm{\Delta }_{m1}^I,$$
(3.22)
for some polynomials $`𝔣(q^\gamma )`$ and $`𝔤(q^\gamma )`$.

We denote by $`\lambda _x`$ and $`\rho _x`$ the operators of left and right multiplication by $`x`$, respectively, on $`\{q^\gamma ,E\}`$, i.e. for $`𝒯\{q^\gamma ,E\}`$ we have

$$\lambda _x(𝒯)=x𝒯\text{and}\rho _x(𝒯)=𝒯x.$$

Eliminating $`q^\gamma `$ (resp. $`q^\gamma `$) from the equations

$`\mathrm{\Lambda }_\gamma E^m`$ $`=(abcdq^{\gamma 1}+q^\gamma 1abcdq^1)E^m,`$
$`E^m\mathrm{\Lambda }_\gamma `$ $`=(q^mabcdq^{\gamma 1}+q^mq^\gamma 1abcdq^1)E^m,`$
one obtains

$$\delta _m^+(E^m)=q^\gamma E^m\text{and}\delta _m^{}(E^m)=q^\gamma E^m,$$
(3.23)
with

$`\delta _m^+`$ $`={\displaystyle \frac{q}{abcd(1q^{2m})}}(\lambda _{\mathrm{\Lambda }_\gamma }q^m\rho _{\mathrm{\Lambda }_\gamma })+\left(1+{\displaystyle \frac{q}{abcd}}\right){\displaystyle \frac{1}{1+q^m}}\mathrm{Id},`$ (3.24a)
and
$`\delta _m^{}`$ $`={\displaystyle \frac{1}{1q^{2m}}}(\lambda _{\mathrm{\Lambda }_\gamma }q^m\rho _{\mathrm{\Lambda }_\gamma })+{\displaystyle \frac{1}{1+q^m}}\left({\displaystyle \frac{abcd}{q}}+1\right)\mathrm{Id}`$
$`={\displaystyle \frac{abcd}{q}}\delta _m^+.`$ (3.24b)
Notice that $`\delta _m^+`$ and $`\delta _m^{}`$ preserve the subalgebra $`\mathrm{\Delta }^I`$. Hence, from (3.22), using the properties (3.23) of $`\delta _m^+`$ and $`\delta _m^{}`$, one deduces that

$$p(\delta _m^+,\delta _m^{})L_{a,b,c,d}^m=(\text{Id}+I)p(q^\gamma ,q^\gamma )\frac{𝔣(q^\gamma )}{𝔤(q^\gamma )}E^m\text{mod}\mathrm{\Delta }_{m1}^I,$$

for any polynomial $`p(q^\gamma ,q^\gamma )[q^\gamma ,q^\gamma ]`$. In particular, if we choose

$$p(q^\gamma ,q^\gamma )=𝔠(q^\gamma )𝔤(q^\gamma )I\left(𝔡(q^\gamma )𝔣(q^\gamma )\right),$$

we get

$$p(\delta _m^+,\delta _m^{})L_{a,b,c,d}^m=𝔡(q^\gamma )𝔣(q^\gamma )I(𝔡(q^\gamma )𝔣(q^\gamma ))(\text{Id}+I)\frac{𝔠(q^\gamma )}{𝔡(q^\gamma )}E^m\text{mod}\mathrm{\Delta }_{m1}^I.$$

Define

$$\mathrm{\Theta }_m=𝔡(q^\gamma )𝔣(q^\gamma )I(𝔡(q^\gamma )𝔣(q^\gamma )).$$

Since $`\mathrm{\Theta }_m`$ is $`I`$-invariant, by Lemma 3.2, $`\mathrm{\Theta }_m𝒦`$. Also, it is clear from the definitions (3.24a) and (3.24b) of $`\delta _m^+`$ and $`\delta _m^{}`$ that $`p(\delta _m^+,\delta _m^{})L_{a,b,c,d}^m`$. Hence, remembering (3.21), we get

$`\mathrm{\Theta }_m^1p(\delta _m^+,\delta _m^{})L_{a,b,c,d}^m`$ $`=(\text{Id}+I){\displaystyle \frac{𝔠(q^\gamma )}{𝔡(q^\gamma )}}E^m\text{mod}\mathrm{\Delta }_{m1}^I`$
$`=T\text{mod}\mathrm{\Delta }_{m1}^I,`$
which establishes (3.20) with $`V_m=p(\delta _m^+,\delta _m^{})L_{a,b,c,d}^m`$. This completes the proof of Theorem 3.3. ∎

## 4. Contiguous relations

From (2.5), we observe that the function $`s_\gamma (aq,b,c,d;a)`$ has a simple product form, since the $`{}_{8}{}^{}W_{7}^{}`$ series defining it reduces to $`1`$. One checks easily that

$$\frac{s_\gamma (aq,b,c,d;a)}{s_{\gamma 1}(aq,b,c,d;a)}=\frac{\phi _\gamma }{\phi _{\gamma 1}}\frac{C_{\gamma ;a,b,c,d}}{A_{\gamma ;a,b,c,d}},$$

with

$$\phi _\gamma =\frac{1}{q^\gamma abcdq^\gamma }.$$
(4.1)
Let us denote by $`P_{a,b,c,d}`$ and $`Q_{a,b,c,d}`$ the difference operators

$`P_{a,b,c,d}`$ $`=\phi _\gamma (E\mathrm{Id}),`$ (4.2)
$`Q_{a,b,c,d}`$ $`=\left(A_{\gamma ;a,b,c,d}\mathrm{Id}C_{\gamma ;a,b,c,d}E^1\right){\displaystyle \frac{1}{\phi _\gamma }}.`$ (4.3)
From (2.1) and (2.2c), it follows that

$$L_{a,b,c,d}\left(a+a^1\right)\mathrm{Id}=Q_{a,b,c,d}P_{a,b,c,d}.$$
(4.4)
Since the kernel of the operator $`Q_{a,b,c,d}`$ is generated by $`s_\gamma (aq,b,c,d;a)`$, it is natural to expect that interchanging the factors $`Q_{a,b,c,d}`$ and $`P_{a,b,c,d}`$ in (4.4), we obtain the operator $`L_{aq,b,c,d}`$. This simply means that $`L_{a,b,c,d}`$ can be obtained as a Darboux transformation from the operator $`L_{aq,b,c,d}`$. The proposition below makes this statement precise.

###### Proposition 4.1.

With the notations above, we have

$$L_{aq,b,c,d}\left(a+a^1\right)\mathrm{Id}=P_{a,b,c,d}Q_{a,b,c,d}.$$
(4.5)
Moreover, up to a factor independent of $`\gamma `$, $`Q_{a,b,c,d}`$ maps $`r_\gamma (aq,b,c,d;z)`$ and $`s_\gamma (aq,b,c,d;z)`$ into $`r_\gamma (a,b,c,d;z)`$ and $`s_\gamma (a,b,c,d;z)`$, respectively. More precisely, the following contiguous relations

$`Q_{a,b,c,d}r_\gamma (aq,b,c,d;z)`$ $`={\displaystyle \frac{1az}{a}}r_\gamma (a,b,c,d;z),`$ (4.6a)
$`Q_{a,b,c,d}s_\gamma (aq,b,c,d;z)`$ $`={\displaystyle \frac{(za)(1az)}{az}}s_\gamma (a,b,c,d;z),`$ (4.6b)
hold. In particular, from (4.6b), it follows that the kernel of $`Q_{a,b,c,d}`$ is spanned by $`s_\gamma (aq,b,c,d;a)`$.

###### Proof.

A straightforward computation shows that

$$a+a^1C_{\gamma +1;a,b,c,d}A_{\gamma ;a,b,c,d}=aq+(aq)^1A_{\gamma ;aq,b,c,d}C_{\gamma ;aq,b,c,d},$$
(4.7)
which combined with (2.1), (2.2), (4.2), and (4.3) gives (4.5). The proof of (4.6) can be easily extracted from the contiguous relations found by Ismail and Rahman in . Let us denote for simplicity by $`𝔯_\gamma (a,b,c,d;z)`$ and $`𝔰_\gamma (a,b,c,d;z)`$ the $`{}_{8}{}^{}W_{7}^{}`$ factor on the right-hand side of formulas (2.4) and (2.5), respectively. Applying \[22, formula (2.3), p. 207\], with $`a`$, $`b`$, $`c`$, $`d`$, $`e`$, $`f`$ replaced by $`bcd/qz`$, $`abcdq^\gamma `$, $`q^\gamma `$, $`b/z`$, $`c/z`$, $`d/z`$, we obtain

$$\begin{array}{cc}& \left(q^\gamma abcdq^{\gamma 1}\right)\left(1\frac{1}{az}\right)𝔯_\gamma (a,b,c,d;z)\hfill \\ & =\left(q^\gamma 1\right)\left(1\frac{bcdq^{\gamma 1}}{z}\right)𝔯_{\gamma 1}(aq,b,c,d;z)\hfill \\ & +\left(1abcdq^{\gamma 1}\right)\left(1\frac{1}{azq^\gamma }\right)𝔯_\gamma (aq,b,c,d;z).\hfill \end{array}$$
(4.8)
From (2.4) and (4.8), one can deduce (4.6a). To get (4.6b), we first substitute $`bcdzq^{2\gamma +1}`$, $`z/a`$, $`bcq^\gamma `$, $`bdq^\gamma `$, $`cdq^\gamma `$, $`q^{\gamma +1}`$ for $`A^2`$, $`A/\lambda `$, $`A/\mu `$, $`A/\nu `$, $`A/\rho `$, $`A/\sigma `$ in \[22, formula (2.17), p. 210\] to obtain

$`𝔰_{\gamma 1}(aq,b,c,d;z)={\displaystyle \frac{(1bcdzq^{2\gamma })(1bcdzq^{2\gamma 1})}{(1bzq^\gamma )(1czq^\gamma )(1dzq^\gamma )(1bcdzq^{\gamma 1})}}\times `$
$`[(1z/a)(1az)𝔰_{\gamma 1}(a,b,c,d;z)`$ (4.9)
$`+{\displaystyle \frac{z}{a}}{\displaystyle \frac{(1abq^\gamma )(1acq^\gamma )(1adq^\gamma )(1abcdq^{\gamma 1})}{(1abcdq^{2\gamma 1})(1abcdq^{2\gamma })}}𝔰_\gamma (aq,b,c,d;z)],`$
which combined with (2.5) gives (4.6b). ∎

###### Remark 4.2.

We can also consider $`L_{aq,b,c,d}`$ as obtained by means of a Darboux transformation from $`L_{a,b,c,d}`$. In this case, the analogous formulas to (4.6a) and (4.6b) are

$`P_{a,b,c,d}r_\gamma (a,b,c,d;z)`$ $`={\displaystyle \frac{az}{z}}r_\gamma (aq,b,c,d;z),`$ (4.10a)
$`P_{a,b,c,d}s_\gamma (a,b,c,d;z)`$ $`=s_\gamma (aq,b,c,d;z).`$ (4.10b)
One way to prove these formulas is to use \[22, formula (2.2), p. 207\] and \[22, formula (2.18), p. 210\], and to proceed as before, adapting the proof of Proposition 4.1. However, they can be easily deduced from Proposition 4.1, by applying $`P_{a,b,c,d}`$ to both sides of (4.6a)-(4.6b), and by using (4.5) and the fact that $`r_\gamma (a,b,c,d;z)`$ and $`s_\gamma (a,b,c,d;z)`$ satisfy (2.9).

## 5. Computing the kernel of the Darboux transform using Sears’ and Watson’s transformation formulas

As we saw in Theorem 3.3, in order to apply the general bispectral Darboux technology, one needs to find an $`I`$-invariant factorization of a constant coefficient polynomial $``$ in $`L_{a,b,c,d}`$. In this section we construct $`I`$-invariant eigenfunctions of $`L_{a,b,c,d}`$, which generate the kernel of the operator $`𝒬`$, in the notations of Theorem 3.1.

###### Proposition 5.1.

Assume that

$$a=dq^\alpha \text{and}z_k=dq^k,k,\alpha ,\mathrm{\hspace{0.33em}0}k\alpha 1.$$
(5.1)
Then,

$`r_\gamma (z_k)`$ $`=(bcq^k,d^2q^\alpha ;q)_k\times `$
$`f_\gamma {}_{4}{}^{}\varphi _{3}^{}[\begin{array}{c}d^2q^k,bcd^2q^{\gamma +\alpha 1},q^\gamma ,q^k\\ bd,cd,d^2q^\alpha \end{array};q,q],`$ (5.2)
$`s_\gamma (z_k)`$ $`=(bc)^{\alpha k1}(d^2q^k,dq^{k+1}/b,dq^{k+1}/c;q)_{\alpha k1}\times `$
$`f_\gamma g_\gamma {}_{4}{}^{}\varphi _{3}^{}[\begin{array}{c}q^{k\alpha +1},bcq^\gamma ,q^{\gamma \alpha +1}/d^2,q^{1\alpha k}/d^2\\ q^{2\alpha }/d^2,bq^{1\alpha }/d,cq^{1\alpha }/d\end{array};q,q],`$ (5.3)
with

$`f_\gamma `$ $`={\displaystyle \frac{q^{\gamma \alpha }}{(bdq^\gamma ,cdq^\gamma ;q)_\alpha }},`$ (5.4)
$`g_\gamma `$ $`=q^{\gamma (\alpha 1)}d^{2\gamma }{\displaystyle \frac{(d^2q^{\gamma +\alpha },bcd^2q^{\gamma +\alpha 1};q)_{\mathrm{}}}{(q^{\gamma +1},bcq^\gamma ;q)_{\mathrm{}}}}.`$ (5.5)
If we make the further assumption that

$$d^2=q^l,l,l1,$$
(5.6)
$`g_\gamma `$ in (5.5) can be rewritten as

$$g_\gamma =\frac{q^{\gamma (l+\alpha 1)}}{(q^{\gamma +1},bcq^\gamma ;q)_{l+\alpha 1}}.$$
(5.7)
###### Proof.

Applying Watson’s transformation formula \[10, (III.18), p. 242\] with $`n`$, $`a`$, $`b`$, $`c`$, $`d`$, $`e`$ replaced by $`k`$, $`bcq^{k1}`$, $`bq^k/d`$, $`cq^k/d`$, $`bcd^2q^{\gamma +\alpha 1}`$, $`q^\gamma `$, respectively, we can rewrite $`r_\gamma (z_k)`$ as

$`r_\gamma (z_k)`$ $`=q^{(\alpha k)\gamma }{\displaystyle \frac{(d^2q^{\gamma +\alpha },bcq^k;q)_k}{(bdq^\gamma ,cdq^\gamma ;q)_\alpha }}{\displaystyle \frac{(q^{\alpha k+1}/d^2;q)_k}{(q^{\gamma \alpha k+1}/d^2;q)_k}}\times `$
$`{}_{4}{}^{}\varphi _{3}^{}[\begin{array}{c}d^2q^k,bcd^2q^{\gamma +\alpha 1},q^\gamma ,q^k\\ cd,bd,d^2q^\alpha \end{array};q,q].`$ (5.8)
Using formula (I.9) on page 233 in , one can easily deduce that

$$\frac{(q^{\alpha k+1}/d^2;q)_k}{(q^{\gamma \alpha k+1}/d^2;q)_k}=q^{\gamma k}\frac{(d^2q^\alpha ;q)_k}{(d^2q^{\gamma +\alpha };q)_k},$$

which shows that (5.8) is equivalent to (5.2). The proof of (5.3) is similar. Let us first replace the parameters $`n`$, $`a`$, $`b`$, $`c`$, $`d`$, $`e`$ in Watson’s formula \[10, (III.18), p. 242\] by $`\alpha k1`$, $`bcd^2q^{2\gamma +k}`$, $`bdq^\gamma `$, $`cdq^\gamma `$, $`bcq^\gamma `$, $`q^{\gamma +1}`$, to obtain the following formula for $`s_\gamma (z_k)`$ in terms of $`{}_{4}{}^{}\varphi _{3}^{}`$ series

$`s_\gamma (z_k)`$ $`=q^{(\alpha +k)\gamma }d^{2\gamma }{\displaystyle \frac{(d^2q^{\gamma +\alpha },bcd^2q^{\gamma +\alpha 1};q)_{\mathrm{}}}{(q^{\gamma +1},bcq^\gamma ;q)_{\mathrm{}}}}{\displaystyle \frac{(d^2q^k;q)_{\alpha k1}}{(bdq^\gamma ,cdq^\gamma ;q)_{k+1}}}\times `$
$`{}_{4}{}^{}\varphi _{3}^{}[\begin{array}{c}q^{k+1},bcq^\gamma ,q^{\gamma +1},q^{\alpha +k+1}\\ cdq^{\gamma +k+1},bdq^{\gamma +k+1},q^{2\alpha }/d^2\end{array};q,q].`$ (5.9)
Now replacing $`n`$, $`a`$, $`b`$, $`c`$, $`d`$, $`f`$ by $`\alpha k1`$, $`bcq^\gamma `$, $`q^{k+1}`$, $`q^{\gamma +1}`$, $`q^{2\alpha }/d^2`$, $`cdq^{\gamma +k+1}`$, $`bdq^{\gamma +k+1}`$ in Sears’ transformation formula (see \[10, (III.15), p. 242\]), we see that (5.9) gives exactly (5.3). ∎

In the case considered above (i.e. $`a=dq^\alpha `$ and $`d^2=q^l`$) formula (3.15) reduces to

$$I(q^\gamma )=\frac{q^{\gamma \alpha +1}}{bcd^2}=\frac{q^{\gamma \alpha l+1}}{bc}.$$
(5.10)
###### Lemma 5.2.

If (5.6) holds, the functions $`r_\gamma (z_k)`$ and $`s_\gamma (z_k)`$, defined by (5.2) and (5.3) are $`I`$-invariant rational functions in $`q^\gamma `$, i.e. we have $`I(r_\gamma (z_k))=r_\gamma (z_k)`$ and $`I(s_\gamma (z_k))=s_\gamma (z_k)`$.

###### Proof.

From (5.2), (5.3), (5.4), and (5.7), it is obvious that $`r_\gamma (z_k)`$ and $`s_\gamma (z_k)`$ are rational functions in $`q^\gamma `$. Since $`I(bcd^2q^{\gamma +\alpha 1})=q^\gamma `$ and $`I(bcq^\gamma )=q^{\gamma \alpha +1}/d^2`$, it is clear that the $`{}_{4}{}^{}\varphi _{3}^{}`$ series in (5.2) and (5.3) are $`I`$-invariant. Thus, to show the invariance of $`r_\gamma (z_k)`$ and $`s_\gamma (z_k)`$, it is enough to show the invariance of $`f_\gamma `$ and $`g_\gamma `$. Using the definitions (5.4) of $`f_\gamma `$ and (5.10) of $`I`$, one computes

$$\begin{array}{cc}\hfill I(f_\gamma )& =\frac{\left(bcd^2q^{\gamma +\alpha 1}\right)^\alpha }{(11/cdq^{\gamma +\alpha 1})(11/cdq^{\gamma +\alpha 2})\mathrm{}(11/cdq^\gamma )}\hfill \\ & \times \frac{1}{(11/bdq^{\gamma +\alpha 1})(11/bdq^{\gamma +\alpha 2})\mathrm{}(11/bdq^\gamma )}\hfill \\ & =\frac{q^{2\gamma \alpha +\alpha (\alpha 1)\alpha (\gamma +\alpha 1)}}{(bdq^\gamma ,cdq^\gamma ;q)_\alpha }=f_\gamma .\hfill \end{array}$$

The proof of the $`I`$-invariance of $`g_\gamma `$ is similar. ∎

## 6. Askey-Wilson type functions

In this section we prove the main result of the paper. First we review some basic facts about Darboux transformations.

### 6.1. Iterating the Darboux transformation

Consider a difference operator $`L_0`$ and a nonzero eigenfunction $`\psi _1(\gamma )`$ for which $`L_0\psi _1(\gamma )=x_1\psi _1(\gamma )`$. If $`Q_0`$ is a first-order difference operator with kernel spanned by $`\psi _1(\gamma )`$, we can write $`L_0`$ in the form $`L_0=x_1\mathrm{Id}+P_0Q_0`$, for some difference operator $`P_0`$. The operator $`L_1=x_1\mathrm{Id}+Q_0P_0`$ is by definition a Darboux transformation of $`L_0`$.

The next proposition, which can be found in \[43, pp. 14-19\] for the case of differential operators, describes the result of $`m`$ successive Darboux transformations starting from $`L_0`$

$`L_0=x_{j_1}\mathrm{Id}+P_0Q_0L_1=x_{j_1}\mathrm{Id}+Q_0P_0=x_{j_2}\mathrm{Id}+P_1Q_1\mathrm{}`$
$`L_{m1}=x_{j_{m1}}\mathrm{Id}+Q_{m2}P_{m2}=x_{j_m}\mathrm{Id}+P_{m1}Q_{m1}`$ (6.1)
$`\widehat{L}=L_m=x_{j_m}\mathrm{Id}+Q_{m1}P_{m1}.`$
###### Proposition 6.1.

If the operator $`\widehat{L}`$ is obtained from $`L_0`$ by iteration of the Darboux transformation (6.1) then

$$\widehat{L}𝒬=𝒬L_0,$$
(6.2)
where

$$𝒬=Q_{m1}Q_{m2}\mathrm{}Q_0.$$
(6.3)
Conversely, if $`\widehat{L}`$ and $`L_0`$ are difference operators and if there exists a difference operator $`𝒬`$ of order $`1`$ such that (6.2) holds, then $`\widehat{L}`$ can be obtained from $`L_0`$ by a sequence of Darboux transformations.

###### Proof.

From (6.1) one can deduce that

$$L_kQ_{k1}=Q_{k1}L_{k1}\text{ for }k=1,2,\mathrm{},m,$$
(6.4)
which easily implies (6.2), thus proving the first assertion.

The second part can be proved by induction on the order of the operator $`𝒬`$ as follows. Assume that (6.2) holds. It is obvious that $`\mathrm{Ker}𝒬`$ is preserved by $`L_0`$, i.e. $`L_0(\mathrm{Ker}𝒬)\mathrm{Ker}𝒬`$. Thus we may regard $`L_0`$ restricted to $`\mathrm{Ker}𝒬`$ as a linear operator acting in a finite-dimensional complex vector space. Let $`\{\psi _j(\gamma )\}_{j=1}^m`$ be a Jordan basis of this operator, i.e.

$$\mathrm{Ker}𝒬=\mathrm{Span}\{\psi _1(\gamma ),\psi _2(\gamma ),\mathrm{},\psi _m(\gamma )\},$$
(6.5a)
where
$$L_0\psi _k(\gamma )=x_{j_k}\psi _k(\gamma )+\sigma _k\psi _{k1}(\gamma ),1km,$$
(6.5b)
with $`\sigma _1=0`$ and $`\sigma _k=0`$ or $`1`$ for $`2km`$. Let $`Q_0`$ be a first-order difference operator with kernel spanned by $`\psi _1(\gamma )`$. If the order of $`𝒬`$ is 1, we simply take $`Q_0=𝒬`$. Otherwise we can take for example

$$Q_0=E\frac{\psi _1(\gamma +1)}{\psi _1(\gamma )}\mathrm{Id}.$$

Since $`\psi _1(\gamma )`$ is in the kernels of $`L_0x_{j_1}\mathrm{Id}`$ and $`𝒬`$, there exist difference operators $`P_0`$ and $`\stackrel{~}{𝒬}`$ such that

$`L_0x_{j_1}\mathrm{Id}=P_0Q_0`$ (6.6)
and
$`𝒬=\stackrel{~}{𝒬}Q_0.`$ (6.7)
If we define $`L_1=x_{j_1}\mathrm{Id}+Q_0P_0`$, then from (6.6) we get that

$$L_1Q_0=Q_0L_0.$$
(6.8)
Using the last equation, (6.2) and (6.7) one can show that

$$\widehat{L}\stackrel{~}{𝒬}=\stackrel{~}{𝒬}L_1.$$

Moreover from (6.5b) and (6.8) we can see that for $`k=2,3,\mathrm{},m`$ we have

$$L_1(Q_0\psi _k)=x_{j_k}Q_0\psi _k+\sigma _k^{}Q_0\psi _{k1},\text{ where }\sigma _2^{}=0,\text{ and }\sigma _k^{}=\sigma _k\text{ if }k>2,$$

i.e. $`\{Q_0\psi _j(\gamma )\}_{j=2}^m`$ is a Jordan basis for the operator $`L_1`$ in the vector space $`\mathrm{Ker}\stackrel{~}{𝒬}`$. The proof now follows by induction. ∎

The above proposition tells us that the Darboux process (6.1) is determined by the intertwining operator $`𝒬`$ in (6.2). The points $`\{x_{j_k}\}`$ are simply the eigenvalues of the operator $`L_0`$ in the vector space $`\mathrm{Ker}𝒬`$. We will use Proposition 6.1 as follows. First, we will construct a difference operator $`𝒬`$ whose kernel is $`L_0`$ invariant, with $`L_0=L_{a,b,c,d}`$ the Askey-Wilson operator defined in Section 2. This operator $`𝒬`$ will correspond to a sequence of Darboux transformations (6.1). Indeed if we consider the operator $`𝒬L_0`$, then its kernel contains the kernel of $`𝒬`$, so it has a factorization of the form

$$\widehat{L}𝒬=𝒬L_0,$$

for some (unique) difference operator $`\widehat{L}`$. By the above considerations, $`\widehat{L}`$ is obtained from $`L_0`$ by a Darboux chain (6.1).

The next important fact is that the operator $`𝒬`$ needed above can be explicitly reconstructed (up to a factor on the right) using the functions $`\{\psi _j\}`$ from its kernel (6.5a). More precisely, if $`[m_1,m_2]`$ is the support of $`𝒬`$, with $`m_2m_1=m`$ we can write the operator $`𝒬`$ in the form

$$𝒬=𝔞(\gamma )\left|\begin{array}{cccc}\psi _1(\gamma +m_2)& \mathrm{}& \psi _m(\gamma +m_2)& E^{m_2}\hfill \\ \mathrm{}& & \mathrm{}& \mathrm{}\hfill \\ \psi _1(\gamma +k)& \mathrm{}& \psi _m(\gamma +k)& E^k\hfill \\ \mathrm{}& & \mathrm{}& \mathrm{}\hfill \\ \psi _1(\gamma +m_1)& \mathrm{}& \psi _m(\gamma +m_1)& E^{m_1}\hfill \end{array}\right|,$$
(6.9)
where $`𝔞(\gamma )`$ is some function in $`\gamma `$, and the determinant is expanded from left to right, i.e. the shift operators $`E^k`$ are pulled to the right. Indeed, it is obvious that the right-hand side of (6.9) defines a difference operator with support $`[m_1,m_2]`$, whose kernel coincides with the kernel of $`𝒬`$. But two operators having the same supports and kernels must differ by a factor, which proves (6.9).

Another useful way to describe the Darboux process (6.1) is to define a new operator $``$ by the following formula

$$=𝔥(L_0)=\underset{i=1}{\overset{m}{}}(L_0x_{j_i}\mathrm{Id}).$$
(6.10)
Using (6.1)-(6.4) one can easily show by induction that

$$=𝒫𝒬,$$
(6.11)
where

$$𝒫=P_0P_1\mathrm{}P_{m1}.$$
(6.12)
Thus, from (6.2) we get

$$\widehat{}:=𝒬𝒫=\underset{i=1}{\overset{m}{}}(\widehat{L}x_{j_i}\mathrm{Id})=𝔥(\widehat{L}),$$
(6.13)
i.e. we have the Darboux map

$$=𝔥(L_0)=𝒫𝒬\widehat{}=𝒬𝒫=𝔥(\widehat{L}).$$
(6.14)
### 6.2. The main result

Now we are ready to formulate our main result. We show in Theorem 6.2 below, that if

$$d=\pm q^{l/2}\text{ and }a=dq^\alpha ,\text{ for some }l,\alpha ,l,\alpha 1,$$
(6.15)
then, the second-order difference operator $`\widehat{L}`$ obtained by iteration of the Darboux process (6.1), starting with the Askey-Wilson operator $`L_0=L_{a,b,c,d}`$, and iterating at any subset $`\{x_{j_1},\mathrm{},x_{j_m}\}`$ of the points

$$x_k=z_k+z_k^1,z_k=dq^k=\pm q^{l/2+k}0k\alpha 1,$$
(6.16)
possesses a two-dimensional space of eigenfunctions which are also eigenfunctions of a $`q`$-difference operator in the spectral variable $`z`$.

Notice that in this case the functions $`\{\psi _j\}`$ must be linear combinations of the functions $`r_\gamma (z_{j_k})`$ and $`s_\gamma (z_{j_k})`$, computed in Proposition 5.1, i.e.

$$\psi _k(\gamma )=\mu _kr_\gamma (a,b,c,d;z_{j_k})+\nu _ks_\gamma (a,b,c,d;z_{j_k}),k=1,2,\mathrm{},m,$$
(6.17)
for some constants $`\{\mu _k,\nu _k\}_{k=1}^m`$. Moreover the operator $`L_0`$ is diagonalizable (in the vector space $`\mathrm{Ker}𝒬`$), i.e. $`\sigma _k=0`$ for all $`k`$ in (6.5b). It is interesting that this is a purely $`q`$-phenomenon. In the $`q=1`$ case, one needs to use generalized eigenvectors, see  and .

###### Theorem 6.2.

Assume that (6.15) holds, i.e. $`a=\pm q^{\frac{l}{2}+\alpha }`$ and $`d=\pm q^{\frac{l}{2}}`$, $`l,\alpha `$, $`l,\alpha 1`$, with the same choice of sign for $`a`$ and $`d`$. The second-order difference operator $`\widehat{L}`$ which is obtained by iteration of the Darboux transformation (6.1), starting with the associated Askey-Wilson operator $`L_0=L_{a,b,c,d}`$ (2.1), and adding bound states at the points

$`\{x_k=dq^k+(dq^k)^1=\pm (q^{\frac{l}{2}+k}+q^{\left(\frac{l}{2}+k\right)}),`$
$`k\{j_1,j_2,\mathrm{},j_m\}\{0,1,\mathrm{},\alpha 1\}\},`$ (6.18)
possesses a two-dimensional space of eigenfunctions which are also eigenfunctions of a $`q`$-difference operator in the spectral variable $`z`$.

If $`\mathrm{\Psi }(\gamma ,z)`$ satisfies (2.9) (i.e. it is a linear combination of $`R_\gamma `$ and $`S_\gamma `$ defined by (2.8)), then the function $`\widehat{\mathrm{\Psi }}(\gamma ,z)`$ defined by

$$\widehat{\mathrm{\Psi }}(\gamma ,z)=𝒬\mathrm{\Psi }(\gamma ,z),$$
(6.19)
satisfies

$$\widehat{L}\widehat{\mathrm{\Psi }}(\gamma ,z)=(z+z^1)\widehat{\mathrm{\Psi }}(\gamma ,z).$$
(6.20)
Equation (6.20) follows immediately from (6.2). So, it remains to show that there exists a $`q`$-difference operator $`\widehat{B}(z,D_z)`$ in $`z`$, such that

$$\widehat{B}(z,D_z)\widehat{\mathrm{\Psi }}(\gamma ,z)=\widehat{\mathrm{\Lambda }}_\gamma \widehat{\mathrm{\Psi }}(\gamma ,z),$$
(6.21)
for some function $`\widehat{\mathrm{\Lambda }}_\gamma `$.

Before we start with the proof of this fact, we shall make two simple but important remarks.

###### Remark 6.3.

So far, we have not said anything about the normalization and the supports of the operators $`\{P_k,Q_k\}_{k=0}^{m1}`$, or equivalently, about the normalization $`𝔞(\gamma )`$ and the support $`[m_1,m_2]`$ of $`𝒬`$ in (6.9). The point is that this is not important for the bispectral relation (6.21). Indeed, if (6.21) holds, and if

$$𝒬^{}=𝔞^{}(\gamma )E^m^{}𝒬\text{ and }\widehat{\mathrm{\Psi }}^{}(\gamma ,z)=𝒬^{}\mathrm{\Psi }(\gamma ,z)=𝔞^{}(\gamma )\widehat{\mathrm{\Psi }}(\gamma +m^{},z),$$
(6.22)
for some function $`𝔞^{}(\gamma )`$ and some $`m^{}`$, then it follows from (6.21) and (6.22) that

$$\widehat{B}(z,D_z)\widehat{\mathrm{\Psi }}^{}(\gamma ,z)=\widehat{\mathrm{\Lambda }}_\gamma ^{}\widehat{\mathrm{\Psi }}^{}(\gamma ,z),$$
(6.23)
where $`\widehat{\mathrm{\Lambda }}_\gamma ^{}=\widehat{\mathrm{\Lambda }}_{\gamma +m^{}}`$, i.e. we still have a $`q`$-difference equation in $`z`$. Thus, the choice of the support $`[m_1,m_2]`$ and of the normalizing function $`𝔞(\gamma )`$ in (6.9) is not essential for the bispectrality.

###### Remark 6.4.

From Proposition 4.1 we know that $`L_{a,b,c,d}`$ is a Darboux transform from $`L_{aq,b,c,d}`$ at $`a+a^1`$. Thus, we can think of $`\widehat{L}`$ also as a Darboux transformation from the operator $`L_{aq,b,c,d}`$ at the points $`\{x_{j_1},x_{j_2},\mathrm{},x_{j_m},a+a^1\}`$. From (4.4), (4.5), and (6.2) one can deduce that

$$\widehat{L}𝒬^{}=𝒬^{}L_{aq,b,c,d}\text{ with }𝒬^{}=𝒬Q_{a,b,c,d}.$$
(6.24)
Let us define functions $`\psi _k^{}(\gamma )`$ for $`k=1,2,\mathrm{},m+1`$ by

$$\psi _k^{}(\gamma )=(1a/z_{j_k})\mu _kr_\gamma (aq,b,c,d;z_{j_k})+\nu _ks_\gamma (aq,b,c,d;z_{j_k})\text{ for }1km,$$
(6.25)
where $`\mu _k`$ and $`\nu _k`$ are the constants in (6.17), and

$$\psi _{m+1}^{}(\gamma )=s_\gamma (aq,b,c,d;a).$$
(6.26)
Using now equations (4.6a)-(4.6b), we see that

$$Q_{a,b,c,d}\psi _k^{}(\gamma )=\frac{(z_{j_k}a)(1az_{j_k})}{az_{j_k}}\psi _k(\gamma )\text{for }1km,$$
(6.27)
which combined with (6.24) gives

$$\mathrm{ker}𝒬^{}=\mathrm{Span}\{\psi _1^{}(\gamma ),\psi _2^{}(\gamma ),\mathrm{},\psi _{m+1}^{}(\gamma )\}.$$
(6.28)
Thus, we can think of $`\widehat{L}`$ as a Darboux transformation from $`L_{aq,b,c,d}`$, and the kernel of the corresponding intertwining operator is given explicitly by (6.25) and (6.26).

###### Proof of Theorem 6.2.

From Remark 6.4, it follows that we can assume $`m`$ is even, i.e. $`m=2\overline{m}`$ for some $`\overline{m}`$. Remark 6.3 tells us that without any restriction, we can fix the support of $`𝒬`$ to be $`[\overline{m},+\overline{m}]`$ (i.e. $`m_2=m_1=\overline{m}`$) in (6.9), and we can choose an appropriate normalization factor $`𝔞(\gamma )`$ in (6.9). Let us define $`𝔞(\gamma )`$ as

$$𝔞(\gamma ):=\{\begin{array}{cc}1\hfill & \text{if }\overline{m}\text{ is even;}\hfill \\ \frac{1}{q^\gamma abcdq^{\gamma 1}}\hfill & \text{if }\overline{m}\text{ is odd.}\hfill \end{array}$$
(6.29)
Notice that

$$I(𝔞(\gamma ))=(1)^{\overline{m}}𝔞(\gamma ),$$
(6.30)
i.e. up to a sign, $`𝔞(\gamma )`$ is an $`I`$-invariant function. From Lemma 5.2, (6.17), and (6.9) it follows that $`𝒬`$ is a difference operator with coefficients, which are rational functions in $`q^\gamma `$, i.e. $`𝒬\{q^\gamma ,E\}`$. From Lemma 5.2 and (6.17) we know that the kernel of $`𝒬`$ consists of $`I`$-invariant functions, i.e.

$$I(\psi _k(\gamma ))=\psi _k(\gamma )\text{ for }1km.$$

From this relation and the definition of the involution $`I`$, one can easily check that

$$I(\psi _k(\gamma +i))=\psi _k(\gamma i)\text{ for any }i.$$
(6.31)
Formulas (6.9), (6.30), (6.31) show that $`𝒬`$ is an $`I`$-invariant operator, i.e. $`𝒬\mathrm{\Delta }^I`$. Since $`L_0=L_{a,b,c,d}\mathrm{\Delta }^I`$, we see immediately from (6.10) that $`\mathrm{\Delta }^I`$, which combined with (6.11) gives that $`𝒫\mathrm{\Delta }^I`$. The proof now follows from Theorem 3.3 and Theorem 3.1. ∎

###### Remark 6.5.

Let us denote

$$h_a(\gamma )=a^\gamma (abq^\gamma ,acq^\gamma ,adq^\gamma ;q)_{\mathrm{}},$$
(6.32)
and define a new second-order difference operator $`𝔏`$ by conjugating the Askey-Wilson operator $`L_{a,b,c,d}`$ in (2.1) with $`h_a(\gamma )`$, i.e.

$$𝔏=h_a(\gamma )^1L_{a,b,c,d}h_a(\gamma ).$$
(6.33)
If we write $`𝔏`$ in the form

$$𝔏=𝔄_\gamma E+𝔅_\gamma \mathrm{Id}+_\gamma E^1,$$
(6.34)
then a straightforward computation shows that the coefficients of $`𝔏`$ are given by the following formulas

$`𝔄_\gamma `$ $`={\displaystyle \frac{h_a(\gamma +1)}{h_a(\gamma )}}A_\gamma ={\displaystyle \frac{1abcdq^{\gamma 1}}{(1abcdq^{2\gamma 1})(1abcdq^{2\gamma })}}`$ (6.35a)
$`_\gamma `$ $`={\displaystyle \frac{h_a(\gamma 1)}{h_a(\gamma )}}C_\gamma =(1q^\gamma )(1abq^{\gamma 1})(1acq^{\gamma 1})(1adq^{\gamma 1})`$
$`\times {\displaystyle \frac{(1bcq^{\gamma 1})(1bdq^{\gamma 1})(1cdq^{\gamma 1})}{(1abcdq^{2\gamma 2})(1abcdq^{2\gamma 1})}}`$ (6.35b)
$`𝔅_\gamma `$ $`=B_\gamma =q^{\gamma 1}{\displaystyle \frac{(1+abcdq^{2\gamma 1})(sq+s^{}abcd)q^{\gamma 1}(1+q)abcd(s+s^{}q)}{(1abcdq^{2\gamma 2})(1abcdq^{2\gamma })}},`$ (6.35c)
where $`s=a+b+c+d`$ and $`s^{}=a^1+b^1+c^1+d^1`$. Notice that the operator $`𝔏`$ is symmetric in $`a,b,c,d`$ (in fact, this was the operator originally introduced by Askey and Wilson, see \[2, p. 5\]). Therefore, we can produce bispectral operators as in Theorem 6.2, if any two of the parameters $`(a,b,c,d)`$ are related in the same manner as $`a`$ and $`d`$ in (6.15).

Theorem 6.2 and Remark 6.5 tell us that if two of the parameters $`a,b,c,d`$ are related by a formula similar to (6.15), we can construct bispectral operators by adding bound states to the left or to the right of the continuous spectrum, depending on the sign of these two parameters (they must have the same signs!). Assume now that we have two positive parameters, satisfying (6.15) with $`+`$, and the other two parameters are negative and satisfy (6.15) with $``$. Then it is possible to construct bispectral operators by adding bound states both to the left and to the right of the continuous spectrum. The corollary below explains what modifications are necessary in the proof of Theorem 6.2 in this case. The operators considered below contain the so called continuous $`q`$-Jacobi operators which can be obtained from the Askey-Wilson’s operator by taking

$$a=q^{\alpha +\frac{1}{2}},b=q^{\beta +\frac{1}{2}},c=q^{\frac{1}{2}},d=q^{\frac{1}{2}}.$$

This substitution was used by Rahman in . For a different parametrization of the continuous $`q`$-Jacobi polynomials, as well as relations between the different parametrizations see \[25, p. 83\].

###### Corollary 6.6.

Assume that

$$a=q^{l_1/2+\alpha },b=q^{l_2/2+\beta },c=q^{l_2/2},d=q^{l_1/2},$$
(6.36)
where $`\alpha ,\beta ,l_1,l_2`$ are positive integers. Let us define

$`J`$ $`=\left\{q^{l_1/2+k}+q^{l_1/2k},\mathrm{\hspace{0.33em}0}k\alpha 1\right\},`$
$`K`$ $`=\left\{\left(q^{l_2/2+k}+q^{l_2/2k}\right),\mathrm{\hspace{0.33em}0}k\beta 1\right\},`$
and let $`X=\{x_{j_1},x_{j_2},\mathrm{},x_{j_m}\}`$ be a subset of $`JK`$. The second-order difference operator $`\widehat{L}`$ which is obtained by iteration of the Darboux transformation (6.1), starting with the associated Askey-Wilson operator $`L_0=L_{a,b,c,d}`$ (2.1), and adding bound states at the points of $`X`$ possesses a two-dimensional space of eigenfunctions which are also eigenfunctions of a $`q`$-difference operator in the spectral variable $`z`$.

###### Proof.

As in the proof of Theorem 6.2, we can assume that $`m`$ is even, i.e. $`m=2\overline{m}`$ and the support of $`𝒬`$ is $`[\overline{m},\overline{m}]`$. From Lemma 5.2 we know that the eigenfunctions of $`L`$, corresponding to eigenvalues $`x_{j_k}J`$ are $`I`$-invariant. The main point is to see what happens with the eigenfunctions $`\psi _k(\gamma )`$, corresponding to $`x_{j_k}K`$. Let us denote by $`𝔏^b`$ the Askey-Wilson operator with parameters $`a,b,c,d`$ replaced by $`b,c,d,a`$, respectively, i.e.

$$𝔏^b=L_{b,c,d,a}.$$
(6.37)
Then using the notations in Remark 6.5 we will have

$$𝔏=h_b(\gamma )^1𝔏^bh_b(\gamma ),$$
(6.38)
where $`h_b(\gamma )`$ is analogous to $`h_a(\gamma )`$ defined in (6.32):

$$h_b(\gamma )=b^\gamma (baq^\gamma ,bcq^\gamma ,bdq^\gamma ;q)_{\mathrm{}}.$$
(6.39)
From (6.33) and (6.38) it follows that

$$L_{a,b,c,d}=\frac{h_a(\gamma )}{h_b(\gamma )}𝔏^b\frac{h_b(\gamma )}{h_a(\gamma )}.$$
(6.40)
This equation simply means that the correspondence between the eigenfunctions of $`𝔏^b`$ and $`L_{a,b,c,d}`$ is given by multiplication by the function $`h_a(\gamma )/h_b(\gamma )`$, i.e. if

$$\mathrm{\Psi }_a(\gamma )=\frac{h_a(\gamma )}{h_b(\gamma )}\mathrm{\Psi }_b(\gamma ),$$

then

$$L_{a,b,c,d}\mathrm{\Psi }_a(\gamma )=\lambda \mathrm{\Psi }_a(\gamma )\text{if and only if}𝔏^b\mathrm{\Psi }_b(\gamma )=\lambda \mathrm{\Psi }_b(\gamma ).$$

If we take an eigenfunction of $`𝔏^b`$ with eigenvalue $`x_{j_k}K`$, then this function must be $`I`$-invariant by Lemma 5.2, interchanging the roles of the parameters $`(a,d)`$ and $`(b,c)`$. Thus, it remains to see how $`h_a(\gamma )/h_b(\gamma )`$ is transformed under the involution $`I`$. Using (6.36), (6.32) and (6.39) we get

$$\frac{h_a(\gamma )}{h_b(\gamma )}=\left(\frac{a}{b}\right)^\gamma \frac{(acq^\gamma ,adq^\gamma ;q)_{\mathrm{}}}{(bcq^\gamma ,bdq^\gamma ;q)_{\mathrm{}}}=(1)^\gamma h(\gamma ),$$
(6.41)
where

$$h(\gamma )=q^{[(l_1l_2)/2+\alpha \beta ]\gamma }\frac{(q^{\gamma +\alpha +l_1},q^{\gamma +(l_1+l_2)/2+\alpha };q)_{\mathrm{}}}{(q^{\gamma +\beta +l_2},q^{\gamma +(l_1+l_2)/2+\beta };q)_{\mathrm{}}}.$$
(6.42)
Using the definitions at the beginning of Section 2 we can write $`h(\gamma )`$ in the form

$$h(\gamma )=\frac{q^{[(l_1l_2)/2+\alpha \beta ]\gamma }}{(q^{\gamma +\beta +l_2};q)_{\alpha \beta +l_1l_2}(q^{\gamma +(l_1+l_2)/2+\beta };q)_{\alpha \beta }},$$
(6.43)
where

$$(y;q)_n=\frac{1}{(yq^n;q)_n}\text{ for }n<0.$$

Notice that if $`l_1l_2`$ is even, then $`h(\gamma )`$ is a rational function in $`q^\gamma `$, i.e. $`h(\gamma )\{q^\gamma ,E\}`$. However, for $`l_1l_2`$ odd, $`h(\gamma )`$ has an extra factor $`q^{\gamma /2}`$, i.e. $`h(\gamma )=q^{\gamma /2}\times `$(rational function in $`q^\gamma `$).

In the case considered here (i.e. assuming (6.36)), the involution $`I`$ acts on rational functions in $`q^\gamma `$ by the following formula

$$I(q^\gamma )=q^{\gamma \alpha \beta l_1l_2+1}.$$

We can naturally extend the involution $`I`$ to an involution acting on $`\{q^{\gamma /2},E\}`$ by defining

$$I(q^{\gamma /2})=q^{(\gamma +\alpha +\beta +l_1+l_21)/2}.$$
(6.44)
Using (6.43) and (6.44) one can check that

$$I(h(\gamma ))=(1)^{\alpha +\beta +l_1+l_2}h(\gamma ),$$
(6.45)
i.e. up to a sign the function $`h(\gamma )`$ is $`I`$-invariant. Finally, let us denote by $`m^{}`$ and $`m^{\prime \prime }`$ the number of elements in the sets $`XJ`$ and $`XK`$, respectively. The computations above can be summarized as follows:

* $`m^{}`$ of the functions $`\psi _k`$, generating the kernel of $`𝒬`$, are eigenfunctions of $`L_{a,b,c,d}`$ corresponding to eigenvalues $`x_{j_k}J`$ and these functions are $`I`$-invariant rational functions in $`q^\gamma `$;
* $`m^{\prime \prime }`$ of the functions $`\psi _k`$, generating the kernel of $`𝒬`$, are eigenfunctions of $`L_{a,b,c,d}`$ corresponding to eigenvalues $`x_{j_k}K`$ and these functions are of the form

$$\psi _k(\gamma )=(1)^\gamma \xi _k(\gamma ),$$

where

$$I(\xi _k(\gamma ))=(1)^{\alpha +\beta +l_1+l_2}\xi _k(\gamma ),$$
(6.46)
and

$$\xi _k(\gamma )=\{\begin{array}{cc}\text{rational function in }q^\gamma \hfill & \text{ if }l_1l_2\text{ is even;}\hfill \\ q^{\gamma /2}\times (\text{rational function in }q^\gamma )\hfill & \text{ if }l_1l_2\text{ is odd.}\hfill \end{array}$$
(6.47)

Without any restriction, we can assume that $`x_{j_k}J`$ for $`k=1,2,\mathrm{},m^{}`$ and $`x_{j_s}K`$ for $`s=m^{}+1,\mathrm{},m^{\prime \prime }`$. Let us denote by $`(a_{ij};b_{ij};c_i)`$ the $`(m+1)\times (m+1)`$ matrix with rows numbered from $`\overline{m}`$ to $`\overline{m}`$ (i.e. $`i=\overline{m},\overline{m}1,\mathrm{},\overline{m}`$) and columns from 1 to $`m+1`$, with entries $`a_{ij}`$ for the first $`m^{}`$ columns, entries $`b_{ij}`$ for the next $`m^{\prime \prime }`$ columns, and entries $`c_i`$ in the last column. Then, from (6.9) we see that the operator $`𝒬`$ has the form

$$𝒬=𝔞(\gamma )(1)^{m^{\prime \prime }\gamma }det(\psi _j(\gamma +i);(1)^i\xi _j(\gamma +i);E^i),$$
(6.48)
where $`𝔞(\gamma )`$ is an appropriate function. From (6.46) it follows that

$$I(\xi _j(\gamma +i))=(1)^{\alpha +\beta +l_1+l_2}\xi _j(\gamma i),$$

which combined with (6.48) shows that

$$I(𝒬)=(1)^MI(𝔞(\gamma )(1)^{m^{\prime \prime }\gamma })det(\psi _j(\gamma +i);(1)^i\xi _j(\gamma +i);E^i),$$
(6.49)
where $`M=\overline{m}+(\alpha +\beta +l_1+l_2)m^{\prime \prime }`$. From (6.47) and (6.48) we see that

$$𝒬=𝔞(\gamma )(1)^{m^{\prime \prime }\gamma }\times \{\begin{array}{cc}\text{operator in }\{q^\gamma ,E\}\hfill & \text{if }m^{\prime \prime }(l_1l_2)\text{ is even;}\hfill \\ q^{\gamma /2}\times (\text{operator in }\{q^\gamma ,E\})\hfill & \text{if }m^{\prime \prime }(l_1l_2)\text{ is odd.}\hfill \end{array}$$
(6.50)
Let us now define $`𝔞(\gamma )`$ as follows

$$𝔞(\gamma )=\frac{1}{(1)^{m^{\prime \prime }\gamma }(q^{ϵ_1\gamma /2}+ϵ_2q^{ϵ_1(\gamma +\alpha +\beta +l_1+l_21)/2})},$$
(6.51)
where

$$ϵ_1=\{\begin{array}{cc}2\hfill & \text{ if }m^{\prime \prime }(l_1l_2)\text{ is even;}\hfill \\ 1\hfill & \text{ if }m^{\prime \prime }(l_1l_2)\text{ is odd;}\hfill \end{array}$$

and

$$ϵ_2=(1)^M=\{\begin{array}{cc}\hfill 1& \text{ if }M\text{ is even;}\hfill \\ \hfill 1& \text{ if }M\text{ is odd.}\hfill \end{array}$$

From (6.48), (6.49), (6.50) it follows that with this definition of $`𝔞(\gamma )`$, $`𝒬`$ becomes an $`I`$-invariant operator from $`\{q^\gamma ,E\}`$, and the proof continues as in Theorem 6.2. ∎

## 7. A $`q`$-analogue of the Krall-Jacobi polynomials, with a bound state off the continuous spectrum

In this section we illustrate Theorem 6.2 on the special case

$$a=q^{\frac{3}{2}},d=q^{\frac{1}{2}},b\text{and}c\text{arbitrary},$$
(7.1)
adding a bound state off the continuous spectrum of $`L_{a,b,c,d}`$ at the point

$$x_0=q^{\frac{1}{2}}+q^{\frac{1}{2}}.$$
(7.2)
This corresponds to the choice

$$l=1,\alpha =1,k=0.$$
(7.3)
### 7.1. Illustrating the general theory

According to Remark 6.4, this case must be handled by first expressing $`L_{a,b,c,d}`$ as a Darboux transform (without free parameter) from $`L_{aq,b,c,d}`$ at $`a+a^1`$, and then performing a Darboux transform (with free parameter) of $`L_{a,b,c,d}`$ at $`x_0`$. Thus we perform the chain of elementary Darboux transformations

$$\begin{array}{c}L_{aq,b,c,d}=(a+a^1)\mathrm{Id}+P_{a,b,c,d}Q_{a,b,c,d}\hfill \\ \hfill L_{a,b,c,d}=(a+a^1)\mathrm{Id}+Q_{a,b,c,d}P_{a,b,c,d}=x_0\mathrm{Id}+PQ\\ \hfill \widehat{L}=x_0\mathrm{Id}+QP.\end{array}$$
(7.4)
We now explain how Theorems 3.1 and 3.3 (which form the core of the general proof of Theorem 6.2) apply to this example. We can factorize

$$=\left(L_{aq,b,c,d}x_0\mathrm{Id}\right)\left(L_{aq,b,c,d}(a+a^1)\mathrm{Id}\right),$$
(7.5)
as

$$=𝒫𝒬\text{with}𝒫=P_{a,b,c,d}P\text{and}𝒬=QQ_{a,b,c,d}.$$
(7.6)
The explicit formulas for $`P_{a,b,c,d}`$ and $`Q_{a,b,c,d}`$ are given in (4.2) and (4.3) respectively. The kernel $`\psi _\gamma `$ of the operator $`L_{a,b,c,d}x_0\mathrm{Id}`$ is readily obtained by substitution of (7.1) and (7.3) in (5.2) and (5.3) (using (5.1), (5.4), (5.6) and (5.7)), which leads to

$$\psi _\gamma =\frac{q^\gamma }{(1bq^{\gamma +\frac{1}{2}})(1cq^{\gamma +\frac{1}{2}})}\left(\mu +\frac{q^\gamma }{(1q^{\gamma +1})(1bcq^\gamma )}\right),$$
(7.7)
with $`\mu `$ a free parameter. This gives

$$Q=\phi _\gamma (\psi _\gamma E\psi _{\gamma +1}),P=\frac{C_{\gamma ;a,b,c,d}}{\psi _\gamma \phi _{\gamma 1}}E^1+\frac{A_{\gamma ;a,b,c,d}}{\psi _\gamma \phi _\gamma }.$$
(7.8)
In this last formula $`\phi _\gamma `$ is as in (4.1) and corresponds to the factor $`𝔞(\gamma )`$ defined in (6.29) (with $`\overline{m}=1`$ and $`a`$ replaced by $`aq`$), and $`A_{\gamma ;a,b,c,d}`$ and $`C_{\gamma ;a,b,c,d}`$ are defined as in (2.2a) and (2.2b); it is understood that the values (7.1) for $`a`$ and $`d`$ are to be substituted everywhere.

The following identities are equivalent to (4.5):

$$\frac{\phi _\gamma }{\phi _{\gamma +1}}A_{\gamma +1;a,b,c,d}=A_{\gamma ;aq,b,c,d}\text{and}\frac{\phi _\gamma }{\phi _{\gamma 1}}C_{\gamma ;a,b,c,d}=C_{\gamma ;aq,b,c,d}.$$

Using these identities, by a straightforward computation, we obtain that

$$\begin{array}{c}𝒬=QQ_{a,b,c,d}=A_{\gamma ;aq,b,c,d}\psi _\gamma E\hfill \\ \hfill (A_{\gamma ;a,b,c,d}\psi _{\gamma +1}+C_{\gamma +1;a,b,c,d}\psi _\gamma )\mathrm{Id}+C_{\gamma ;aq,b,c,d}\psi _{\gamma +1}E^1,\end{array}$$
(7.9)
and

$$\begin{array}{c}𝒫=P_{a,b,c,d}P=\frac{A_{\gamma ;aq,b,c,d}}{\psi _{\gamma +1}}E\hfill \\ \hfill \left(\frac{A_{\gamma ;a,b,c,d}}{\psi _\gamma }+\frac{C_{\gamma +1;a,b,c,d}}{\psi _{\gamma +1}}\right)\mathrm{Id}+\frac{C_{\gamma ;aq,b,c,d}}{\psi _\gamma }E^1.\end{array}$$
(7.10)
Since $``$ in (7.5) is a constant coefficient polynomial in $`L_{aq,b,c,d}`$, the relevant involution $`I`$ is obtained by replacing $`a`$ with $`aq`$ in (3.15), that is

$$I(q^\gamma )=\frac{1}{abcdq^\gamma }=\frac{1}{bcq^{\gamma +2}}\text{and}I(E)=E^1.$$

We know from (3.18) that $`I(A_{\gamma ;aq,b,c,d})=C_{\gamma ;aq,b,c,d}`$. Using (2.2a), (2.2b) and (7.7), it is straightforward to check that

$$I(A_{\gamma ;a,b,c,d})=C_{\gamma +1;a,b,c,d}\text{and}I(\psi _\gamma )=\psi _{\gamma +1}.$$
(7.11)
This makes clear that the operators $`𝒬`$ and $`𝒫`$ as expressed in (7.9) and (7.10) are $`I`$-invariant. To alleviate the notations, we shall write $`L`$ and $`\mathrm{\Lambda }`$ to denote the following difference operators:

$$LL_{aq,b,c,d},\mathrm{\Lambda }q^\gamma (1q^\gamma )(1abcdq^\gamma )\mathrm{Id},\text{with}a,b,c,d\text{as in}(\text{7.1}).$$

According to Theorem 3.3, we must have that

$$𝒫=U\mathrm{\Gamma }^1\text{and}𝒬=\mathrm{\Theta }^1V,$$

with $`U,V=L,\mathrm{\Lambda }`$ and $`\mathrm{\Gamma },\mathrm{\Theta }𝒦=\mathrm{\Lambda }`$.

Using the symbolic software MACSYMA, we have obtained the following explicit formulas for $`\mathrm{\Gamma },\mathrm{\Theta },U`$ and $`V`$. It follows from (7.11) that the function $`\psi _\gamma \psi _{\gamma +1}`$ is $`I`$-invariant, hence it can be written as a rational function of $`\mathrm{\Lambda }`$. Explicitly:

$$\psi _\gamma \psi _{\gamma +1}\mathrm{Id}=\mathrm{\Gamma }(\mathrm{\Lambda })\mathrm{\Theta }^1(\mathrm{\Lambda }),\text{with}$$

$$\mathrm{\Gamma }(\mathrm{\Lambda })=\left[(q1)(bc1)\mu +1\right]\left[(q^21)(bcq1)\mu +q\right]\mathrm{Id}$$

$$+\mu \left[q+1+\mu (q1)(2bcqq+bc2)\right]\mathrm{\Lambda }+\mu ^2\mathrm{\Lambda }^2,$$

$$\mathrm{\Theta }(\mathrm{\Lambda })=q^1\left[\mathrm{\Lambda }+(q^21)(bc1)\mathrm{Id}\right]\left[\mathrm{\Lambda }+(q1)(bcq1)\mathrm{Id}\right]$$

$$\times \left[\mathrm{\Lambda }+(bq^{1/2}1)(cq^{3/2}1)\mathrm{Id}\right]\left[\mathrm{\Lambda }+(bq^{3/2}1)(cq^{1/2}1)\mathrm{Id}\right].$$

Then, one finds that

$`U=`$ $`(r_1\mathrm{\Lambda }+r_2\mathrm{\Lambda }^2+r_3\mathrm{\Lambda }^3)L+L(s_1\mathrm{\Lambda }+s_2\mathrm{\Lambda }^2+s_3\mathrm{\Lambda }^3)`$
$`+tL+u_0\mathrm{Id}+u_1\mathrm{\Lambda }+u_2\mathrm{\Lambda }^2+u_3\mathrm{\Lambda }^3,`$
$`V=`$ $`(v_1\mathrm{\Lambda }+v_2\mathrm{\Lambda }^2+v_3\mathrm{\Lambda }^3)L+L(w_1\mathrm{\Lambda }+w_2\mathrm{\Lambda }^2+w_3\mathrm{\Lambda }^3)`$
$`+x\mathrm{\Lambda }L\mathrm{\Lambda }+y_0\mathrm{Id}+y_1\mathrm{\Lambda }+y_2\mathrm{\Lambda }^2+y_3\mathrm{\Lambda }^3,`$
for some appropriate choice of the constants $`r_i,s_i,t,u_i,v_i,w_i,x,y_i`$. The explicit expressions for these constants are rather complicated and there is no point in displaying them. Theorem 3.1 (see especially (3.10) and (3.11)) provides the following explicit bispectral operator

$$\widehat{B}=𝔟(V)𝔟(U)f^1,$$
(7.12)
with

$`f=`$ $`𝔟()=\left(z+z^1(q^{\frac{1}{2}}+q^{\frac{1}{2}})\right)\left(z+z^1(q^{\frac{3}{2}}+q^{\frac{3}{2}})\right),`$
$`𝔟(U)=`$ $`(z+z^1)(r_1B+r_2B^2+r_3B^3)+(s_1B+s_2B^2+s_3B^3)(z+z^1)`$
$`+t(z+z^1)+u_0+u_1B+u_2B^2+u_3B^3,`$
$`𝔟(V)=`$ $`(z+z^1)(v_1B+v_2B^2+v_3B^3)+(w_1B+w_2B^2+w_3B^3)(z+z^1)`$
$`+xB(z+z^1)B+y_0+y_1B+y_2B^2+y_3B^3,`$
where $`B=B_{q^{5/2},b,c,q^{1/2}}(z,D_z)`$ and $`B_{a,b,c,d}(z,D_z)`$ is defined as in (2.6). Since $`B`$ is a second-order $`q`$-difference operator, $`\widehat{B}`$ in (7.12) is a $`q`$-difference operator of order 12.

### 7.2. A bispectral operator of order 4

It is possible for the example discussed above, with $`a,b,c,d`$ as in (7.1), to produce a bispectral operator of order 4, by performing a bispectral Darboux transformation of the operator

$$=\left(L_{q^{1/2},b,c,q^{1/2}}(q^{1/2}+q^{1/2})\mathrm{Id}\right)^2,$$

instead of the choice we made in (7.5). The idea is to first obtain the operator $`L_{q^{3/2},b,c,q^{1/2}}`$ as a (contiguous) Darboux transform of $`L_{q^{1/2},b,c,q^{1/2}}`$ at $`q^{1/2}+q^{1/2}`$ using (4.4) and (4.5), and then to perform a (general) upper-lower Darboux transform (with a free parameter) of $`L_{q^{3/2},b,c,q^{1/2}}`$ at $`q^{1/2}+q^{1/2}`$. In the case $`q=1`$, such a strategy has already been explained with all details in Section 4 of . Since we don’t know yet how to produce a bispectral operator of the lowest possible order for the general situation described in Theorem 6.2, we shall only display the final result. Let us define

$$D_2(z)=\frac{\mu q^{3/2}(zb)(zc)(zbq)(zcq)(zq^{5/2})(q^{1/2}z1)}{(z+q^{1/2})(z+q^{3/2})(z^21)(z^2q^2)},$$

$$D_1(z)=q(q+1)(zb)(zc)(zq^{1/2})$$

$$\times \frac{(z+q^{3/2})(q^{1/2}z+1)\mu q^{1/2}\left[(q+bc)(z^2+q)(b+c)(q^2+1)z\right]}{(z+q^{1/2})(z+q^{3/2})(q^{1/2}z+1)(z^21)},$$

$$D_1(z)=D_1(1/z),D_2(z)=D_2(1/z),$$

$$D_0(z)=D_2(z)D_1(z)D_1(z)D_2(z),$$
(7.13)
and

$$\widehat{\mathrm{\Lambda }}(\gamma )=\frac{(q^{\gamma +1}1)(bcq^{\gamma +1}1)\left[q^\gamma (q+1)+\mu (bcq^\gamma 1)(q^{\gamma +2}1)\right]}{q^{2\gamma }}.$$

Then, the functions

$$\widehat{\mathrm{\Psi }}(\gamma ,z)=Q\mathrm{\Psi }(\gamma ,z),\text{with}$$

$$\mathrm{\Psi }(\gamma ,z)=R_\gamma (q^{3/2},b,c,q^{1/2};z)\text{or}S_\gamma (q^{3/2},b,c,q^{1/2};z),$$

and with $`Q`$ as in (7.8), $`R_\gamma (a,b,c,d;z)`$ and $`S_\gamma (a,b,c,d;z)`$ as in (2.8a) and (2.8b), satisfy the pair of equations

$$\widehat{L}\widehat{\mathrm{\Psi }}(\gamma ,z)=(z+1/z)\widehat{\mathrm{\Psi }}(\gamma ,z),\widehat{L}\text{as in}(\text{7.4}),$$

$$\underset{i=2}{\overset{2}{}}D_i(z)\widehat{\mathrm{\Psi }}(\gamma ,q^iz)=\widehat{\mathrm{\Lambda }}(\gamma )\widehat{\mathrm{\Psi }}(\gamma ,z).$$

Let us consider the special case of (7.1) given by

$$a=q^{3/2},d=q^{1/2},b=q^{1/2},c=q^{\beta +1/2},\beta >1,$$
(7.14)
and let us also specialize $`\gamma `$ to be a positive integer. Then, the functions

$$R_n(q^{3/2},q^{1/2},q^{\beta +1/2},q^{1/2};z),n=0,1,2,\mathrm{},$$

become (up to a normalizing factor) the continuous $`q`$-Jacobi polynomials with parameters $`(\alpha =1,\beta )`$ as defined by Rahman (see , formula (3.10.14)). When $`q1`$, these polynomials reduce to the standard Jacobi polynomials $`P_n^{(1,\beta )}(x),x=(z+z^1)/2`$ (see , formula (5.10.2)), with weight function $`(1x)(1+x)^\beta `$ on $`[1,1]`$. The new functions

$$\widehat{P}_0(z)=1,$$

$$\widehat{P}_n(z)(QR)_{n1}=\varphi _{n1}\psi _{n1}R_n(q^{3/2},q^{1/2},q^{\beta +1/2},q^{1/2};z)$$

$$\varphi _{n1}\psi _nR_{n1}(q^{3/2},q^{1/2},q^{\beta +1/2},q^{1/2};z),n=1,2,3,\mathrm{},$$

obtained from the Darboux process, with $`Q`$ as in (7.8), are (Laurent) polynomials of degree $`n`$ in the variable $`z+z^1`$. These polynomials are eigenfunctions of the fourth-order $`q`$-difference operator with coefficients as in (7.13), with $`b`$ and $`c`$ as in (7.14).

Putting the free parameter $`\mu `$ of the Darboux transform in (7.7) to be

$$\mu =\frac{\rho }{2^{\beta +1}(q1)^2},$$

one can check that, taking the limit $`q1`$, the polynomials

$$K_n(x)=\underset{q1}{lim}(q1)^3\widehat{P}_n(z),x=\frac{z+z^1}{2},$$

are the so-called Krall-Jacobi polynomials, which are orthogonal on the interval $`[1,1]`$ for the weight function

$$(1+x)^\beta +\rho \delta (x1),$$
(7.15)
where $`\delta (x)`$ denotes Dirac’s delta function. This is one of the family of orthogonal polynomials satisfying a fourth-order differential equation that was discovered by H.L. Krall in . Thus, the example discussed in this section can be viewed as a $`q`$-deformation of the Krall-Jacobi polynomials, with a bound state of arbitrary weight $`\rho `$ off the continuous spectrum, precisely at the point $`x(z+z^1)/2=(q^{1/2}+q^{1/2})/2`$. When $`q1`$, this bound state tends to the boundary $`x=1`$ of the continuous spectrum of the Jacobi polynomials, in agreement with (7.15).

## 8. Askey-Wilson type solitons

In this section we construct a large family of second-order difference operators which provide solutions to the Askey-Wilson bispectral problem as defined in (1.1) and (1.2), within the context of the theory of commutative rank 1 rings of difference operators. We start from the simple observation that when $`a=b=1`$ and $`c=d=\sqrt{q}`$, the Askey-Wilson operator (2.1) reduces to

$$L_{a,b,c,d}=E+E^1.$$
(8.1)
More precisely, we shall prove the next theorem:

###### Theorem 8.1.

Let $`k_i,1ig`$, denote $`g`$ positive integers. Let us partition these integers into two arbitrary sets $`J`$ and $`K`$ of $`\alpha `$ and $`\beta `$ elements respectively, $`\alpha ,\beta 0,\alpha +\beta =g`$. We build accordingly the $`\alpha `$ points $`q^{k_i/2}+q^{k_i/2},iJ`$, and the $`\beta `$ other points $`(q^{k_i/2}+q^{k_i/2}),iK`$, and we assume that all these points are distinct.

i) The tridiagonal operator $`\widehat{L}`$ which is obtained by iteration of the Darboux process starting at $`L=E+E^1`$, and adding bound states at these $`g=\alpha +\beta `$ distinct points, is part of a commutative rank 1 ring of difference operators, with spectral curve as in (1.6).

ii) The (reduced) wave function $`\mathrm{\Psi }(n,z)`$ of this rank 1 commutative ring of difference operators, is also a common eigenfunction of a rank 1 commutative ring of $`q`$-difference operators in the spectral variable z.

The proof of this result is closely connected with the idea of ”duality” that we exploited to establish Theorem 2.1. In the context of rank 1 commutative rings of differential operators, ”duality” was first systematically exploited by G. Wilson . Before proving Theorem 8.1, we need to give a brief summary of the theory of rank 1 commutative rings of difference operators. We limit ourselves to the case when the spectrum of the ring is a rational irreducible curve, which is the only part of the theory that we shall use. The presentation follows closely our previous work .

### 8.1. Commutative rings of difference operators and rational curves

We denote by $`𝒫=[z]`$ the space of polynomials in the variable $`z`$, and we denote by $`=(z)`$ the space of rational functions of $`z`$.

###### Definition 8.2.

A flag $`𝒱:\mathrm{}V_{n+1}V_nV_{n1}\mathrm{}`$ of subspaces $`V_n`$ is called a rational flag if and only if there are polynomials $`p(z)`$ and $`r(z)`$ (independent of n) such that

$$z^np(z)𝒫V_nz^nr^1(z)𝒫,$$

and the codimension of $`V_n`$ in $`z^nr^1(z)𝒫`$ is equal to the degree of $`r(z)`$.

There is a one-to-one correspondence between rational flags and affine irreducible rational curves which complete by adding *two* non-singular points at infinity. The curve is recovered from the flag as the spectrum of the ring

$$\begin{array}{c}A_𝒱=\{f(z)\text{with poles only at}z=0\text{and}z=\mathrm{},\hfill \\ \hfill \text{such that}k\text{for which}f(z)V_nV_{n+k},n\}.\end{array}$$
(8.2)
When the curve is rationally parametrized by $`z`$, the two non-singular distinguished points at infinity $`P_{\mathrm{}}`$ and $`Q_{\mathrm{}}`$ correspond respectively to $`z=\mathrm{}`$ and $`z=0`$.

By definition, the *rank* of a commutative ring of difference operators is the greatest common divisor of the orders of all the operators in the ring. There is a dictionary between the affine rings of irreducible algebraic curves which complete by adding *two* non-singular points at infinity and the rank 1 commutative rings of difference operators. We refer the reader to  for a full account of the theory, as well as for adequate references. We limit ourselves to the case of irreducible rational curves as described above. The general case proceeds along similar lines. To each $`V_n`$ there corresponds a tau function $`\tau _{V_n}(t)\tau (n,t)`$ in the sense of Segal and Wilson , depending on infinitely many variables $`t=(t_1,t_2,t_3,\mathrm{})`$. The *wave* and *adjoint wave functions*, denoted respectively by $`w(n,t,z)`$ and $`w^{}(n,t,z)`$, are defined by

$`w(n,t,z)`$ $`=z^n\text{exp}(t,z){\displaystyle \frac{\tau (n,t[z^1])}{\tau (n,t)}},`$ (8.3)
$`w^{}(n,t,z)`$ $`=z^n\text{exp}^1(t,z){\displaystyle \frac{\tau (n,t+[z^1])}{\tau (n,t)}},\text{ with}`$ (8.4)
$`\mathrm{exp}(t,z)`$ $`=\mathrm{exp}(t_1z+t_2z^2+t_3z^3+\mathrm{})\text{ and }[z]=(z,z^2/2,z^3/3,\mathrm{}).`$
The plane $`V_n`$ can be recovered from the wave function as

$$V_n=\text{Span}\{w(n,0,z),w(n+1,0,z),w(n+2,0,z),\mathrm{}\}.$$

The *dual flag* of subspaces $`𝒱^{}:\mathrm{}V_{n+1}^{}V_n^{}V_{n1}^{}\mathrm{}`$, with

$$V_n^{}=\mathrm{Span}\{w^{}(n,0,z),w^{}(n1,0,z),w^{}(n2,0,z),\mathrm{}\},$$

defines the same tau function, up to a change of sign, $`\tau _{V_n^{}}(t)=\tau _{V_n}(t)`$.

###### Theorem 8.3.

*(see ,  and references therein)*. The wave function $`w(n,t,z)`$ is the common eigenfunction of a commutative rank 1 ring of difference operators which is isomorphic to $`A_𝒱`$, as defined in (8.2). More precisely, for any $`fA_𝒱`$, there is a difference operator $`L_f`$ such that

$$L_fw(n,t,z)=f(z)w(n,t,z).$$

If $`f`$ has a pole of order $`i`$ at $`P_{\mathrm{}}`$ and a pole of order $`j`$ at $`Q_{\mathrm{}}`$, the operator $`L_f`$, thought of as a finite band matrix, has $`i`$ diagonals above the main diagonal and $`j`$ diagonals below it.

We illustrate the concepts above on the example of the simplest rational singular curves, namely those which are obtained by identifying $`2g`$ distinct points on the Riemann sphere $`^1()`$ in pairs $`\{\lambda _i,\mu _i\},i=1,\mathrm{},g`$. We assume that all these points are distinct from the two distinguished points $`Q_{\mathrm{}}=0`$ and $`P_{\mathrm{}}=\mathrm{}`$. We pick $`g`$ nonzero arbitrary complex numbers $`\delta _i^{},1ig`$ (i.e. a divisor of degree $`g`$ on the curve), and we define

$$\begin{array}{c}V_n=\frac{1}{_{i=1}^g(z\lambda _i)}\{\text{meromorphic functions }f\text{ on }^1()\text{ such that}\hfill \\ \hfill (f)nQ_{\mathrm{}}0\text{ on }^1()\{\mathrm{}\}\text{ and }f(\lambda _i)=\delta _if(\mu _i),i=1,\mathrm{},g\}.\end{array}$$
(8.5)
By the notation $`(f)nQ_{\mathrm{}}0`$, we mean that $`f`$ has a zero at least of order $`n`$ at $`Q_{\mathrm{}}`$ if $`n0`$ and a pole of order at most $`n`$ at $`Q_{\mathrm{}}`$ if $`n<0`$. The next lemma shows that this definition fits within Definition 8.2. In order to establish it, it is useful to introduce the following definitions. A function $`f`$ which admits an expansion around $`z=\mathrm{}`$ of the form

$$f=c_nz^n+c_{n1}z^{n1}+\mathrm{},c_n0,$$
(8.6)
will be called an *element of order n*. It will also be useful to consider the non-degenerate bilinear form on $``$

$$B(f,g)=\text{res}_{z=\mathrm{}}f(z)g(z)\text{d}z,f,g.$$
(8.7)
###### Lemma 8.4.

The plane $`V_n`$ defined in (8.5) satisfies

$$z^n\underset{i=1}{\overset{g}{}}(z\mu _i)𝒫V_nz^n\underset{i=1}{\overset{g}{}}(z\lambda _i)^1𝒫,$$

and the codimension of $`V_n`$ in $`z^n_{i=1}^g(z\lambda _i)^1𝒫`$ is equal to $`g`$.

###### Proof.

The second inclusion as well as the assertion about the codimension both follow immediately from the definition of $`V_n`$ in (8.5). In order to establish the first inclusion, we introduce the dual flag of subspaces

$$\begin{array}{c}V_n^{}=\frac{1}{_{i=1}^g(z\mu _i)}\{\text{meromorphic functions }f\text{ on }^1()\text{ such that}\hfill \\ \hfill (f)+nQ_{\mathrm{}}0\text{ on }^1()\{\mathrm{}\}\text{ and }f(\mu _i)=\delta _i^{}f(\lambda _i),i=1,\mathrm{},g\},\end{array}$$
(8.8)
with

$$\delta _i^{}=\delta _i\underset{\begin{array}{c}j=1\\ ji\end{array}}{\overset{g}{}}\frac{(\mu _i\mu _j)(\mu _i\lambda _j)}{(\lambda _i\mu _j)(\lambda _i\lambda _j)}.$$
(8.9)
It follows immediately from (8.5), (8.8) and (8.9) that

$$\{\text{res}_{z=\lambda _i}+\text{res}_{z=\mu _i}\}w(z)w^{}(z)\text{d}z=0,w(z)V_n\text{and}w^{}(z)V_m^{}.$$

For $`w(z)V_n`$ and $`w^{}(z)V_m^{}`$ with $`nm`$, $`w(z)w^{}(z)`$ has no pole at $`z=0`$. Thus, by the residue theorem, we deduce that

$$\begin{array}{c}\text{res}_{z=\mathrm{}}w(z)w^{}(z)\text{d}z=\underset{i=1}{\overset{g}{}}\{\text{res}_{z=\lambda _i}+\text{res}_{z=\mu _i}\}w(z)w^{}(z)\text{d}z=0,\hfill \\ \hfill w(z)V_n\text{and}w^{}(z)V_m^{}\text{with}nm,\end{array}$$

showing that $`V_n^{}V_n^{}`$, with $`V_n^{}`$ the orthogonal of $`V_n`$ with respect to the bilinear form $`B`$ defined in (8.7). Since $`V_n`$ has a basis of elements of orders $`n,n+1,n+2,\mathrm{}`$ (see (8.6) for the definition), and $`V_n^{}`$ has a basis of elements of orders $`n,n+1,n+2,\mathrm{}`$, we deduce that $`V_n^{}=V_n^{}`$. By the definition (8.8) of $`V_n^{}`$, it implies that

$$V_n^{}=V_n^{}z^n\underset{i=1}{\overset{g}{}}(z\mu _i)^1𝒫V_nz^n\underset{i=1}{\overset{g}{}}(z\mu _i)𝒫,$$

as desired. This finishes the proof of the lemma. ∎

###### Lemma 8.5.

The tau function $`\tau (n,t)`$ associated with the rational flag $`𝒱`$ defined by (8.5), is given by

$$\begin{array}{c}\tau (n,t)=\underset{j=1}{\overset{g}{}}\mu _j^n\text{exp}(\underset{i=1}{\overset{\mathrm{}}{}}t_i\underset{j=1}{\overset{g}{}}\lambda _j^i)\times \hfill \\ \hfill \text{det}\left(\lambda _j^{n+i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{n+i1}\text{exp}(t,\mu _j)\right)_{1i,jg}.\end{array}$$
(8.10)
###### Proof.

It follows easily from the definition of $`V_n`$ in (8.5) that the wave function associated with this plane is

$$\begin{array}{c}w(n,t,z)=\frac{\text{exp}(t,z)}{_{j=1}^g(z\lambda _j)}\times \hfill \\ \hfill \frac{\text{det}(\lambda _j^{n+i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{n+i1}\text{exp}(t,\mu _j);z^{n+i1})}{\text{det}\left(\lambda _j^{n+i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{n+i1}\text{exp}(t,\mu _j)\right)_{1i,jg}}.\end{array}$$
(8.11)
The notation in the numerator means the determinant of the $`(g+1)\times (g+1)`$ matrix with $`(i,j)`$ entries as indicated there for $`1ig+1,1jg`$, and entries $`z^{n+i1},1ig+1`$, in the last column. Elementary row manipulations with this determinant lead to

$$w(n,t,z)=\frac{\text{exp}(t,z)}{_{j=1}^g(z\lambda _j)}z^{n+g}\frac{_{j=1}^g(z\lambda _j)\tau (n,t[z^1])}{z^g\tau (n,t)},$$

in agreement with (8.3), which defines the tau function. This finishes the proof. ∎

###### Remark 8.6.

The motivation for the (irrelevant) factor $`_{j=1}^g\mu _j^n`$ in the definition (8.10) of $`\tau (n,t)`$ will appear in the next subsection, see (8.36).

Having in view Theorem 8.1, we are particularly interested in rational curves with double points, for which the associated commutative rings of difference operators contain a tridiagonal operator.

###### Proposition 8.7.

Let us assume that $`\mu _i=\lambda _i^1`$, and that the 2g points $`\lambda _i,\mu _i,1ig`$, are still distinct of each other, i.e. $`\lambda _i\pm 1`$ and $`\lambda _i\lambda _j^{\pm 1},ij`$. Then, the rank 1 commutative ring of difference operators $`A_𝒱`$, with $`V_n`$ as in (8.5), contains a tridiagonal operator $`\widehat{L}`$ with one diagonal above and one diagonal below the main diagonal. The spectral curve of the ring has for equation

$$\text{Spec}(A_𝒱):v^2=(u^24)\underset{i=1}{\overset{g}{}}(u(\lambda _i+\lambda _i^1))^2,$$
(8.12)
and the operator $`\widehat{L}`$ can be obtained by iterating the Darboux process starting with $`L=E+E^1`$ as in (8.1), and adding bound states at the $`g`$ points $`\lambda _i+\lambda _i^1`$.

###### Proof.

Since $`\mu _i=\lambda _i^1`$, it is clear that the function

$$u=z+z^1,$$
(8.13)
belongs to the ring $`A_𝒱`$ as defined in (8.2), with $`uV_nV_{n1}`$. This function has a simple pole both at $`z=0`$ and at $`z=\mathrm{}`$ hence, by Theorem 8.3, there exists a tridiagonal operator $`L_u`$ with one diagonal above and one diagonal below the main diagonal, satisfying $`L_uw(n,t,z)=uw(n,t,z)`$.

Since $`u`$ has a simple pole at $`z=0`$, all other generators of $`A_𝒱`$ can be taken to be polynomials. The functions

$$f_k=z^k\underset{i=1}{\overset{g}{}}(z\lambda _i)(z\lambda _i^1),k,$$
(8.14)
which vanish at $`\lambda _i`$ and $`\mu _i=\lambda _i^1`$, obviously belong to $`A_𝒱`$. For $`1kg1`$, by subtracting an appropriate polynomial of degree $`k`$ in $`u`$, one obtains polynomials (in $`z`$) $`q_{2g1}(z),q_{2g2}(z),\mathrm{},q_{g+1}(z)`$, of degrees $`2g1,2g2,\mathrm{},g+1`$, that belong to $`A_𝒱`$. Combining them with the functions $`f_k,k0`$, defined in (8.14), we conclude that $`A_𝒱`$ contains polynomials $`q_k(z)`$ of degree $`k`$, $`kg+1`$. On the other hand, any polynomial $`q(z)A_𝒱`$ of degree $`kg`$, must be identically constant. Indeed, such a polynomial must satisfy $`q(\lambda _i)=q(\lambda _i^1),\mathrm{\hspace{0.33em}1}ig`$, which amounts to a linear homogeneous system $`Ac=0`$, for the unknown coefficients $`c_1,\mathrm{},c_g`$ of $`q(z)=c_0+c_1z+\mathrm{}+c_gz^g`$. One shows easily that

$$\text{det}A=\underset{i=1}{\overset{g}{}}\lambda _i^g\underset{i=1}{\overset{g}{}}(\lambda _i^21)\underset{1i<jg}{}(\lambda _i\lambda _j)\underset{1i<jg}{}(1\lambda _i\lambda _j).$$

Since the $`2g`$ points $`\lambda _i,\mu _i=\lambda _i^1,1ig`$, are assumed to be distinct, this determinant is nonzero, implying that $`c_1=c_2=\mathrm{}=c_g=0`$.

The upshot of the discussion above is that the algebra $`A_𝒱`$ is generated by $`u`$ in (8.13) and polynomials of degree $`k,kg+1`$, or equivalently by $`u`$ and the functions $`f_k,k`$, introduced in (8.14). We now show that the functions $`u`$ and

$$v=\frac{zz^1}{z^g}\underset{i=1}{\overset{g}{}}(z\lambda _i)(z\lambda _i^1),$$

are enough to generate $`A_𝒱`$, which will establish (8.12). Clearly

$$f_g=\underset{i=1}{\overset{g}{}}\left(u(\lambda _i+\lambda _i^1)\right)\text{and}v=f_{g1}f_{g+1}.$$

Since $`uf_k=f_{k1}+f_{k+1}`$, we deduce inductively that $`f_{g1},f_{g+1},f_{g2},f_{g+2},\mathrm{}`$, belong to the algebra generated by $`u`$ and $`v`$.

It is well known that the spectral meaning of the curve (8.12) is that the operator $`L_u\widehat{L}`$, $`u`$ as in (8.13), can be obtained by iteration of the Darboux process, starting from $`L=E+E^1`$ and adding bound states at the points $`\lambda _i+\lambda _i^1`$. This can be easily deduced from the explicit form for the wave function given in (8.11), with $`\mu _i=\lambda _i^1`$. The proof is complete. ∎

### 8.2. The proof of Theorem 8.1

The proof of Theorem 8.1 is a combination of Proposition 8.7 above and Proposition 8.8 below. It is convenient to introduce the so-called reduced wave and adjoint wave functions, which are obtained by omitting the factors $`\text{exp}(t,z)`$ and $`\text{exp}^1(t,z)`$ in (8.3) and (8.4) respectively:

$$\mathrm{\Psi }(n,z)=\text{exp}^1(t,z)w(n,t,z),\mathrm{\Psi }^{}(n,z)=\text{exp}(t,z)w^{}(n,t,z).$$
(8.15)
We also omit to write the explicit dependence on $`t`$, which is irrelevant for what follows; we just think of $`t_1,t_2,t_3,\mathrm{}`$ as some free parameters.

###### Proposition 8.8.

Consider a rational curve with double points, obtained by identifying $`2g`$ distinct points of $`^1()`$ in pairs $`\{\lambda _i,\mu _i\},i=1,\mathrm{},g`$. Assume that

$$\lambda _i=q^{k_i}\mu _i,\text{with}k_i\{0\}.$$
(8.16)
Then the reduced wave function $`\mathrm{\Psi }(n,z)`$, besides being a common eigenfunction of a rank 1 commutative ring $`A_𝒱`$ of difference operators (in $`n`$)

$$L_f\mathrm{\Psi }(n,z)=f(z)\mathrm{\Psi }(n,z),fA_𝒱,$$

is also a common eigenfunction of a rank 1 commutative ring of $`q`$-difference operators $`A_{\stackrel{~}{𝒱}}`$, in the variable z

$$B_g\mathrm{\Psi }(n,z)=g(n)\mathrm{\Psi }(n,z),gA_{\stackrel{~}{𝒱}}.$$
(8.17)
Assuming for a moment this result we can prove Theorem 8.1:

###### Proof of Theorem 8.1.

From (8.16) we immediately obtain that the only intersection between Proposition 8.7 and Proposition 8.8 is when

$$\lambda _i=q^{k_i}\mu _i\text{and}\mu _i=\lambda _i^1\lambda _i=\pm q^{k_i/2},k_i\{0\},$$

in which case (8.12) reduces to (1.6). ∎

As already announced, the proof of Proposition 8.8 exploits an interesting duality of the (reduced) wave function.

###### Proof of Proposition 8.8.

Using (8.16), we have that

$$\begin{array}{c}\lambda _j^{n+i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{n+i1}\text{exp}(t,\mu _j)=\hfill \\ \hfill \mu _j^n\left\{x^{k_j}\lambda _j^{i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{i1}\text{exp}(t,\mu _j)\right\},\text{with}x=q^n.\end{array}$$
(8.18)
From (8.11) it follows that the (reduced) wave function is a function of $`x=q^n`$ and $`z`$, i.e.

$$\begin{array}{c}\mathrm{\Psi }(n,x)\mathrm{\Psi }(x,z)=\hfill \\ \hfill \frac{e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}}}{_{j=1}^g(z\lambda _j)}\frac{\text{det}(x^{k_j}\lambda _j^{i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{i1}\text{exp}(t,\mu _j);z^{i1})}{\text{det}\left(x^{k_j}\lambda _j^{i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{i1}\text{exp}(t,\mu _j)\right)_{1i,jg}}.\end{array}$$
(8.19)
Any (reduced) wave function admits an expansion

$$\mathrm{\Psi }(x,z)=e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}}\left(1+\frac{w_1(x)}{z}+\frac{w_2(x)}{z^2}+\mathrm{}\right)\text{as}z\mathrm{}.$$
(8.20)
Since the $`k_i`$ are positive integers, it follows immediately from (8.19) that

$$\underset{x\mathrm{}}{lim}e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}}\mathrm{\Psi }(x,z)=\frac{\text{det}(\lambda _j^{i1};z^{i1})}{_{j=1}^g(z\lambda _j)\text{det}(\lambda _j^{i1})_{1i,jg}}=1,$$
(8.21)
showing that, when the special conditions (8.16) are satisfied, we also have

$$\mathrm{\Psi }(x,z)=e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}}\left(1+\frac{\stackrel{~}{w}_1(z)}{x}+\frac{\stackrel{~}{w}_2(z)}{x^2}+\mathrm{}\right)\text{as}x\mathrm{}.$$

This suggests that the function $`\stackrel{~}{\mathrm{\Psi }}(x,z)`$ obtained by exchanging the variables $`x`$ and $`z`$ in $`\mathrm{\Psi }(x,z)`$

$$\stackrel{~}{\mathrm{\Psi }}(x,z)=\mathrm{\Psi }(z,x),$$
(8.22)
should be the (reduced) wave function associated with *another* rational flag $`\stackrel{~}{𝒱}:\mathrm{}\stackrel{~}{V}_{n+1}\stackrel{~}{V}_n\stackrel{~}{V}_{n1}\mathrm{}`$, in the sense of Definition 8.2. The correctness of this assertion will be established in Lemma 8.9 below. Assuming the result, we deduce from Theorem 8.3 that

$$B_g\stackrel{~}{\mathrm{\Psi }}(x,z)\underset{\text{finitely many }i\text{’s}}{}D_i(x)\stackrel{~}{\mathrm{\Psi }}(q^ix,z)=g(z)\stackrel{~}{\mathrm{\Psi }}(x,z),g(z)A_{\stackrel{~}{𝒱}},$$

which, remembering the definitions of $`\mathrm{\Psi }(x,z)`$ and $`\stackrel{~}{\mathrm{\Psi }}(x,z)`$ in (8.19) and (8.22), amounts to

$$\underset{\text{finitely many }i\text{’s }}{}D_i(z)\mathrm{\Psi }(n,q^iz)=g(q^n)\mathrm{\Psi }(n,z),$$

which establishes (8.17). This finishes the proof of Proposition 8.8. ∎

As mentioned in the proof of Proposition 8.8 given above, to complete the argument, we still need to establish that $`\stackrel{~}{\mathrm{\Psi }}(x,z)`$ defined as in (8.22) is indeed the (reduced) wave function of some rational flag. We need some preliminaries.

We introduce the following multiplicative groups of formal pseudodifference operators

$$𝒲=\left\{1+\underset{i,j=1}{\overset{\mathrm{}}{}}w_{ij}x^jD_x^i\right\},𝒲^{}=\left\{1+\underset{i,j=1}{\overset{\mathrm{}}{}}w_{ij}x^jD_x^i\right\},$$
(8.23)
with $`w_{ij}`$, $`D_x`$ and $`D_x^1`$ the forward and backward $`q`$-shift operators, $`D_xh(x)=h(qx)`$ and $`D_x^1h(x)=h(q^1x)`$.

We denote by $`𝔞`$ the adjoint isomorphism

$$𝔞:𝒲𝒲^{}:W𝔞(W)=(W^1)^{},$$
(8.24)
where  denotes the adjoint operator, $`(D_x^i)^{}=D_x^i`$. We define an anti-isomorphism $`𝔟:𝒲𝒲`$ by

$$\begin{array}{c}𝔟(x)=D_x,𝔟(D_x)=x,\text{i.e.}\hfill \\ \hfill 𝔟\left(1+\underset{i,j=1}{\overset{\mathrm{}}{}}w_{ij}x^jD_x^i\right)=1+\underset{i,j=1}{\overset{\mathrm{}}{}}w_{ij}x^iD_x^j,\end{array}$$
(8.25)
as well as an anti-isomorphism $`𝔟^{}:𝒲^{}𝒲^{}`$ by

$$\begin{array}{c}𝔟^{}(x)=D_x^1,𝔟^{}(D_x)=x^1,\text{i.e.}\hfill \\ \hfill 𝔟^{}\left(1+\underset{i,j=1}{\overset{\mathrm{}}{}}w_{ij}x^jD_x^i\right)=1+\underset{i,j=1}{\overset{\mathrm{}}{}}w_{ij}x^iD_x^j.\end{array}$$
(8.26)
With these definitions, it is straightforward to check that

$$𝔟^{}𝔞=𝔞𝔟.$$
(8.27)
We can write (8.20) as

$$\begin{array}{c}\mathrm{\Psi }(x,z)=W(x)e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}},\hfill \\ \hfill W(x)=1+w_1(x)D_x^1+w_2(x)D_x^2+\mathrm{}.\end{array}$$
(8.28)
The functions $`w_i(x)`$ are rational functions of $`x`$ which, remembering (8.21), satisfy $`lim_x\mathrm{}w_i(x)=0`$. Thus

$$W(x)=1+\underset{i,j=1}{\overset{\mathrm{}}{}}w_{ij}x^jD_x^i\text{as}x\mathrm{},$$
(8.29)
which shows that $`W(x)𝒲`$ as defined in (8.23).

It is a general fact from the theory of the discrete KP-hierarchy (see , ) that the adjoint (reduced) wave function as defined by (8.4) and (8.15) can be expressed in terms of $`W(x)`$ as follows

$$\mathrm{\Psi }^{}(x,z)=𝔞(W)(q^1x)e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}},$$
(8.30)
with $`𝔞(W)`$ defined as in (8.24). In the same papers, one can find a proof that (8.28) and (8.30) imply the so-called *bilinear identities*

$$B(\mathrm{\Psi }(q^ix,z),\mathrm{\Psi }^{}(x,z))=0,i0,$$
(8.31)
where $`B`$ is the residue pairing introduced in (8.7).

With these preliminaries, we can establish the last needed lemma.

###### Lemma 8.9.

With the assumptions of Proposition 8.8, the new functions

$$\stackrel{~}{\mathrm{\Psi }}(x,z)=\mathrm{\Psi }(z,x)\text{and}\stackrel{~}{\mathrm{\Psi }}^{}(x,z)=\mathrm{\Psi }^{}(qz,q^1x)q^1xz^1,$$
(8.32)
are respectively the (reduced) wave and adjoint wave functions associated with a rational flag $`\stackrel{~}{𝒱}`$.

###### Proof.

In terms of the anti-isomorphism $`𝔟`$ introduced in (8.25), $`\stackrel{~}{\mathrm{\Psi }}(x,z)`$ as defined in (8.32) can be written as follows

$$\stackrel{~}{\mathrm{\Psi }}(x,z)=𝔟(W)(x)e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}},$$
(8.33)
with $`W(x)`$ as in (8.29). Using (8.30) and the definitions (8.26) and (8.32) of $`𝔟^{}`$ and $`\stackrel{~}{\mathrm{\Psi }}^{}(x,z)`$ respectively, one checks easily that

$$\stackrel{~}{\mathrm{\Psi }}^{}(x,z)=𝔟^{}\left(𝔞(W)\right)(q^1x)e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}},$$

from which, using (8.27), we deduce that

$$\stackrel{~}{\mathrm{\Psi }}^{}(x,z)=𝔞\left(𝔟(W)\right)(q^1x)e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}}.$$
(8.34)
Let us define

$`\stackrel{~}{V}_n`$ $`=\text{Span}\{\stackrel{~}{\mathrm{\Psi }}(n,z),\stackrel{~}{\mathrm{\Psi }}(n+1,z),\stackrel{~}{\mathrm{\Psi }}(n+2,z),\mathrm{}\}\text{and}`$
$`\stackrel{~}{V}_n^{}`$ $`=\text{Span}\{\stackrel{~}{\mathrm{\Psi }}^{}(n,z),\stackrel{~}{\mathrm{\Psi }}^{}(n1,z),\stackrel{~}{\mathrm{\Psi }}^{}(n2,z),\mathrm{}\},`$
with $`\stackrel{~}{\mathrm{\Psi }}(n,z)`$ and $`\stackrel{~}{\mathrm{\Psi }}^{}(n,z)`$ the functions obtained by substituting $`q^n`$ for $`x`$ in $`\stackrel{~}{\mathrm{\Psi }}(x,z)`$ and $`\stackrel{~}{\mathrm{\Psi }}^{}(x,z)`$ respectively. As recalled above (see (8.28), (8.30) and (8.31), with $`W`$ replaced by $`𝔟(W)`$), the equations (8.33) and (8.34) imply that $`\stackrel{~}{\mathrm{\Psi }}(n,z)`$ and $`\stackrel{~}{\mathrm{\Psi }}^{}(n,z)`$ satisfy the bilinear identities $`B(\stackrel{~}{\mathrm{\Psi }}(n,z),\stackrel{~}{\mathrm{\Psi }}^{}(m,z))=0,nm`$. Thus

$$\stackrel{~}{V}_n^{}=\stackrel{~}{V}_n^{},$$
(8.35)
where $`\stackrel{~}{V}_n^{}`$ denotes the orthogonal of $`\stackrel{~}{V}_n`$ with respect to $`B`$.

Using (8.18), we can write $`\tau (n,t)`$ in (8.10) as

$$\begin{array}{c}\tau (n,t)\tau (x,t)=\text{exp}(\underset{i=1}{\overset{\mathrm{}}{}}t_i\underset{j=1}{\overset{g}{}}\lambda _j^i)\times \hfill \\ \hfill \text{det}\left(x^{k_j}\lambda _j^{i1}\text{exp}(t,\lambda _j)\delta _j\mu _j^{i1}\text{exp}(t,\mu _j)\right)_{1i,jg},\end{array}$$
(8.36)
with $`x=q^n`$, and thus $`\tau (x,t)`$ is a polynomial in $`x`$. Equations (8.3) and (8.4) combined with the definition (8.32) of $`\stackrel{~}{\mathrm{\Psi }}(x,z)`$ and $`\stackrel{~}{\mathrm{\Psi }}^{}(x,z)`$, show then that

$$\stackrel{~}{V}_nz^n\tau (z,t)^1𝒫\text{and}\stackrel{~}{V}_n^{}z^n\tau (qz,t)^1𝒫,$$

which, because of (8.35), is equivalent to

$$z^n\tau (qz,t)𝒫\stackrel{~}{V}_nz^n\tau (z,t)^1𝒫.$$

Clearly, $`\stackrel{~}{V}_n`$ contains elements of orders $`n,n+1,n+2,\mathrm{}`$, and only those, meaning that the codimension of $`\stackrel{~}{V}_n`$ in $`z^n\tau (z,t)^1𝒫`$ is equal to the degree (as a polynomial in $`z`$) of $`\tau (z,t)`$. Thus the flag $`\stackrel{~}{𝒱}:\mathrm{}\stackrel{~}{V}_{n+1}\stackrel{~}{V}_n\stackrel{~}{V}_{n1}\mathrm{}`$ is a rational flag in the sense of Definition 8.2, and $`\stackrel{~}{\mathrm{\Psi }}(n,z)`$ and $`\stackrel{~}{\mathrm{\Psi }}^{}(n,z)`$ are the (reduced) wave and adjoint wave functions associated with this flag. The lemma is established. ∎

###### Remark 8.10.

Proposition 8.8 has been obtained by F.W. Nijhoff and O.A. Chalykh  in the special case when $`\lambda _i=q\mu _i`$, by writing the (reduced) wave function as

$$\mathrm{\Psi }(x,z)=e^{\frac{\mathrm{log}x\mathrm{log}z}{\mathrm{log}q}}\text{det}\left\{\mathrm{Id}(x\mathrm{Id}+X)^1[X,Z]_{q^1}(z\mathrm{Id}+Z)^1\right\},$$

with $`X`$ and $`Z`$ some $`g\times g`$ matrices, so that the $`q`$-commutator $`[X,Z]_{q^1}`$ is a rank 1 matrix. This argument is insufficient to deal with the (non-generic) situation that concerns us in Theorem 8.1. They conjecture that all bispectral rank 1 commutative rings of difference operators can be parametrized in this way. In the limit $`q=1`$, corresponding to rank 1 bispectral rings of differential operators, this was proved in a highly non-trivial paper by G. Wilson . P. Iliev  has obtained independently of  similar formulas.