{'fraction_non_alphanumeric': 0.0673076923076923, 'fraction_numerical': 0.022584973166368514, 'mean_word_length': 3.425108225108225, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 20, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 22, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Discrete orthogonal polynomials and difference equations of several variables

The goal of this work is to characterize all second order difference operators of several variables that have discrete orthogonal polynomials as eigenfunctions. Under some mild assumptions, we give a complete solution of the problem.


Introduction
The goal of this study is to characterize all second order difference equations of several variables that have discrete orthogonal polynomials as eigenfunctions. More precisely, we consider difference operators of the form where △ i and ∇ i are the forward and backward operator in the direction of the i-th coordinate of R d , respectively, I is the identity operator, A i,j , B i and C are functions of x ∈ R d . The discrete orthogonal polynomials in our consideration are polynomials that are orthogonal with respect to an inner product of the form where V is a lattice set in R d and W is some positive weight function on V . There is a close correlation between D, W and V . Some restrictions need to be imposed on V due to the complexity of the geometry in higher dimensions. Under some mild and, we believe, reasonable assumptions on V , we give a complete solution of the problem.
For d = 1, the one dimensional case, the classification problem was studied by several authors early in the last century, we refer to [1,6] for references. It was found that the classical discrete orthogonal polynomials, namely, Hahn polynomials, Meixner polynomials, Krawtchouk polynomials, and Charlier polynomials are eigenfunctions of second order difference operators on the real line, and these are believed to be the only ones that have such property. To our great surprise, however, another family of solutions turns up when we analyze the problem carefully. The difference equation satisfied by this family of solutions is similar to that of Hahn polynomials but the parameters need to be chosen in a different way. In other words, the difference equation has two separate families of solutions when the parameters are chosen differently.
The same phenomenon also appears in the case of two variables. The second order difference equations that have orthogonal polynomials as eigenfunctions were identified in [12], but not all families of orthogonal polynomials were found. In fact, viewing it as an analogue of second order differential operators that have orthogonal polynomials as eigenfunctions, only one family of solutions that had been known in the literature was identified for each difference equation. The solutions of difference equations, however, turn out to be far richer than that of differential equations. For example, in the case of quadratic eigenvalues, only Hahn polynomials on the set V = {(x, y) : x ≥ 0, y ≥ 0, x + y ≤ N } are identified in [12]; there are in fact several other families, including orthogonal polynomials on N 2 0 , on V = {(x, y) : 0 ≤ x ≤ N 1 , 0 ≤ y ≤ N 2 }, and several others. This will be shown in the present paper.
We start with identifying difference operators D that are self-adjoint with respect to the inner product ·, · . This leads to compatibility conditions between the coefficients of the operator D. The requirement that the difference equations have polynomial solutions and some mild restrictions on V reduce the coefficients of D to special simple forms, which allows us to use the compatibility conditions to determine the coefficients of the difference operator. A careful analysis of the result determines all possible solutions. It is well known that orthogonal bases in several variables are not unique [2]. For each solution, we give a family of mutually orthogonal polynomials explicitly and compute their norm; in other words, we give a family of orthonormal basis explicitly.
It should be mentioned that the analysis in several variables is by no means a straightforward extension of analysis in one or two variables. The complexity of the problem increases substantially as the dimension grows. For example, the number of solutions grows exponentially with the dimension, and the geometry of the admissible lattice sets becomes more involved as the dimension increases. Although our approach resembles the one used in [12], several new ideas are introduced to resolve various outstanding problems. The study in this paper is also much more systematic and thorough, as demonstrated by the new results obtained even in the case of d = 2 and by the new family of discrete orthogonal polynomials of one variable.
The paper is organized as follows. The self-adjoint operators are treated in Section 2, where the compatibility conditions are derived. In Section 3, we study difference equations that have polynomial solutions, called admissible equations, and show that the eigenvalues of such equations are either a quadratic polynomial or a linear polynomial of the index, and the coefficients of such equations need to have certain simple forms under some mild restrictions on V . The compatibility conditions can then be used to determine the coefficients, thus D, explicitly. The case of one variable is discussed in Section 4, in which the new family of discrete orthogonal polynomials is treated in details. For several variables, the case of quadratic eigenvalues is studied in Section 5 and the case of linear eigenvalues is developed in Section 6.

Self-adjoint difference operators
In this section we discuss the basic properties of self-adjoint second order difference operators with respect to an inner product defined on a discrete set. We start with some notations.
Let V be a subset of R d and let W : V → R be a positive function. Consider the inner product defined by the weight W as follows for functions f, g : R d → R. We will be mainly interested in the case when f (x) and g(x) are polynomials of x ∈ R d , but the general discussion in this section does not depend essentially on the class of functions. Denote by {e 1 , e 2 , . . . , e d } the standard basis for R d . We denote by E i , △ i and ∇ i , respectively, the customary shift, forward and backward difference operators acting on a function f (x) as follows where I is the identity operator. We will consider second order difference operators of the form where A i,j , B i and C are some functions of x. Immediately from the definition one can see that D can also be rewritten as where the new coefficients α i,j , β i , γ i and δ are related to the old ones via the formulas The representation (2.2) will be more convenient when we deal with polynomials, because the operators △ i and ∇ i decrease the total degree of a polynomial by one. However, necessary and sufficient conditions for an operator D to be self-adjoint with respect to the inner product (2.1) are much simpler and natural if we write it as in (2.3). Formulas (2.4)-(2.7) allow us to go easily from one representation to another. Let us define "directional" boundaries of V as follows: Example 2.1. Let N 0 denote the set of all nonnegative integers. For the set Other sets that play crucial role later are listed below where l i are positive integers; where N is a positive integer; where S is a nonempty subset of {1, 2, . . . , d} and l i are integers such that 1 and Writing similar equalities for , and comparing Du, v with u, Dv gives the stated conditions.
As an immediate corollary from the relations established above we can deduce compatibility conditions that need to be satisfied by the coefficients of the operator D.
Corollary 2.3 (Compatibility conditions). If the operator D is self-adjoint with respect to the inner product defined by (2.1) then for 1 ≤ i = j ≤ d we have and using (2.10) on the left side and (2.9) for the right side we obtain (2.12a). Similarly if we use (2.9) for all ratios in the identity we get (2.12b), which completes the proof.
Remark 2.4. The compatibility conditions stated above can be easily extended to more general difference operators. For example, we can consider the operator D that adds two additional terms E i E j and E −1 i E −1 j in the operator D in (2.3). However, the complexity of the computations in the sections below will increase significantly even with these two terms added. Furthermore, since △ = ∇ + ∇△, the operator D includes the 4-th order term △ 1 ∇ 1 △ 2 ∇ 2 (say, d = 2) which is the discrete analog of ∂ 2 1 ∂ 2 2 , where ∂ i stands for the partial derivative with respect to x i . Classifying orthogonal polynomials that are eigenfunctions of D appears to be an interesting but much harder problem.

Admissible equations and orthogonal polynomials
In this section we study difference operators D that have discrete orthogonal polynomials as eigenfunctions. For more details on the general theory of discrete orthogonal polynomials of several variables we refer the reader to [11].
Throughout the paper we use the standard multi-index notation. A multi-index will be denoted by µ = (µ 1 , µ 2 , . . . , µ d ) ∈ N d 0 . For each µ we denote by x µ the monomial The degree of a polynomial is defined as the highest degree of its monomials. We denote by R[x] = R[x 1 , x 2 , . . . , x d ] the space of all polynomials in the variables x 1 , x 2 , . . . , x d and by Π d n the subspace of polynomials of degree at most n in the variables {x 1 , x 2 , . . . , x d }. The latter has dimension dim(Π d n ) = n+d n . Let V be an at most countable set of isolated points in R d and let |V | denote the cardinality of V . If f, g = 0 for the inner product defined by (2.1), we say that f and g are orthogonal with respect to W on the set V . Orthogonal polynomials on V depend on the structure of the ideal In fact, the orthogonal polynomials belong to the space There is a lattice set Λ(V ) such that every polynomial in R[V ] can be written as The lattice set is not uniquely determined by V . However, in all cases that arise here, there is a natural way to identify V and Λ(V ), see Example 3.1. Denote Λ k (V ) = {µ ∈ Λ(V ) : |µ| = k} and let r k be the number of the elements in the set Λ k (V ).
. For all other sets defined in Example 2.1 one can easily find (inductively on d and |V |) explicit sets of generators of the ideal I(V ). In particular, the generators listed below allow us to identify Λ(V ) and V .
For V d l , defined by (2.8a), the ideal I(V d l ) is generated by the set , the corresponding ideal is generated by Next we want to consider difference operators D on the space R[V ]. In order to have a proper action of D on R[V ] the ideal I(V ) must be D-invariant, i.e.
One can easily show that if D is self-adjoint with respect to the inner product (2.1), then the above condition is satisfied. More precisely, we have the following proposition. Proof. Notice that if a polynomial p(x) vanishes on V , then the polynomial p(x−e j ) will vanish on V \ ∂ − j V . Thus, if p(x) ∈ I(V ) and (2.11a) holds then γ j (x)p(x − e j ) ∈ I(V ), i.e. the operator γ j E −1 j preserves I(V ). Similar arguments show that (2.11b) (resp. (2.11c)) implies that β j E j (resp. α i,j E i E −1 j ) preserves I(V ), which completes the proof. is called admissible on V if for any k ∈ N 0 there is a number λ k such that the equation Du = λ k u has r k linearly independent polynomial solutions in R[V ] and it has no nontrivial solutions in the set of polynomials of degree less than k.
Remark 3.4. Notice that if (3.1) is admissible on V and if we have nontrivial solutions of the equations Du = λ k u and Du = λ l u for distinct integers k, l ∈ N 0 then λ k = λ l . Thus if D is also a self-adjoint operator with respect to the inner product (2.1), the polynomial solutions for Du = λ k u and Du = λ l u will be mutually orthogonal.
From now on, we will assume that λ 0 = 0 (which is equivalent to C = 0 for operators D of the form (2.2)). This is not a real restriction on D since one can replace D by D − λ 0 I if necessary.
The next proposition shows that we can pick linearly independent polynomial solutions in such a way that the highest terms contain single monomials.
The proof of Proposition 3.5 follows easily form Definition 3.3, see [12] in the case d = 2.
Remark 3.6. If D is an admissible operator on V we can trivially extend it to an admissible operator on the set When we classify the possible admissible operators, we will exclude such degenerate situations.
The next proposition essentially characterizes the admissible difference operators of the form (2.2). In order to eliminate several singular cases we will assume that µ ∈ Λ(V ) for |µ| = 1, 2, 3. Instead of putting µ ∈ Λ(V ) for |µ| = 3 we can require that 3e i ∈ Λ(V ) for some i ∈ {1, 2, . . . , d} and for all 1 ≤ i = j ≤ d at least one of {2e i + e j , e i + 2e j } belongs to Λ(V ), but this would make the statement of the theorem long and awkward. Theorem 3.7. Assume that µ ∈ Λ(V ) for |µ| = 1, 2, 3. Let D be a second order difference operator defined by (2.2) such that I(V ) is D-invariant. Then the following conditions are equivalent: (i) The equation Du = λu is admissible; (ii) The coefficients B i (x) are polynomials of degree at most 1, and A i,j (x) are polynomials of degree at most 2 satisfying for some constants a and b. The eigenvalues λ k must be distinct and they are given by Proof. Assume first that the equation Du = λu is admissible. Applying Proposition 3.5 we see that there exist polynomials Similarly, we can find polynomials of the form , which, combined with (3.6) and (3.7), shows that Next, for i = j we take a polynomial of the form . A simple computation as above gives that shows that (3.3c) holds. It remains to show that A i,j are polynomials of degree at most 2 and that (3.4) holds. First we use a polynomial of the form for some i satisfying Du = λ 3 u, (3.8) to show that (3.4) holds for k = 3. Indeed, we already know that if P is a polynomial of degree at most 2, D(P ) will also be a polynomial of degree at most 2 (because we can write it as a linear combination of 1, P ei , P ei+ej , for i, j = 1, 2, . . . , d). Thus where in the last equality we used (3.3a) and (3.3b). Comparing now the coefficients of x 3 i on both sides in (3.8), we see that λ 3 = 6a + 3b, which is exactly (3.4) for k = 3.
Next we consider the solution to equation (3.8) of the form which combined with (3.3c) proves that A i,j and A j,i are polynomials of degree at most 2. Finally, let µ ∈ Λ(V ) such that |µ| = k and let The admissibility implies that and therefore which combined with (3.9) gives (3.4) Conversely, assume that (ii) holds. One can easily show by induction on |µ| that for µ ∈ N d there exist polynomials of the form satisfying (3.9) with eigenvalue λ |µ| given by (3.4). Indeed, let us assume that this is true for |µ| ≤ k − 1 and take µ such that |µ| = k. Using computation (3.11) we see that for some constants γ µ,ν . Then i.e. if we pick c µ,ν = −γ µ,ν /(λ |ν| − λ |µ| ) the polynomial P µ defined by (3.12) will satisfy (3.9). Thus, in the quotient space R[V ] we have r k polynomials of total degree k satisfying (3.9).
In the next sections we will use conditions (3.3) and the compatibility conditions (2.11) and (2.12) to determine the possible self-adjoint operators on appropriate sets V .
Since the difference equation Du = λu is invariant under translations x → x + h, we can consider the difference equations modulo translations.
In the rest of the paper we will impose certain boundary conditions on the difference operator D. More precisely, we will assume that for every j the coefficients A i,j (x) vanish on the hyperplane {x : x j = 0}. Roughly speaking, this means that V has sufficiently many boundary points belonging to {x : x j = 0}. We show below, that if a = 0 and if V has "enough" points on its boundary then, up to a translation, this must be true.
Let us denote by We say that V has a nontrivial boundary if for every j = 1, 2, . . . , d, ∂ j V contains more than 2 d points and that these points do not belong to a variety of dimension less than d − 1. For example, a line in d = 2 is a nontrivial boundary, so is a plane in d = 3.
Proposition 3.8. Let V be a discrete set with a nontrivial boundary and let D be an admissible, self-adjoint operator. If a = 0, then after an appropriate translation, ∂ j V becomes a subset of {x : x j = 0} and the coefficients A i,j (x) vanish on {x : Proof. Using the fact that D is self-adjoint, Proposition 2.2, (2.4) and (2.6) we see that Applying Bezout's theorem for polynomials in d variables [7, Chapt. IV, Section 2] for the polynomials A 1,j (x), . . . , A d,j (x) which satisfy (3.3b)-(3.3c) and vanish on ∂ j V , one concludes that they must contain a common linear factor. Equation (3.3b) shows that this factor must be of the form (3.13) and completes the proof.
The compatibility conditions also imply that A i,j − A j,i ∈ Π d 1 for i = j, which combined with (3.3c) shows that A i,j = ax i x j mod Π d 1 (the argument works simultaneously for a = 0 and a = 0). The form of the coefficients determined below will be the staring point in the next sections. Proposition 3.9. Let A i,j and B i satisfy the conditions in Theorem 3.7(ii) and assume that A i,j vanishes on {x : x j = 0}. Then, the polynomial identities (2.12a) imply that (3.14) or equivalently: for some constants a, b, l i,j and s i .
Proof. Equation (3.15) follows immediately from (3.3a). From equations (3.3b), (3.3c) and the fact that A i,j vanishes on {x : x j = 0} we deduce that Comparing the coefficient of x 3 i x 3 j on both sides of (2.12a) we obtain In the rest of the paper, we will find the general solution of the compatibility conditions (2.12), by equating the coefficients of the different powers of x, and we will determine the possible sets V , weights W as well as explicit bases of orthogonal polynomials.
Notice that (3.4) implies that λ k must be at most quadratic in k. Thus, we have essentially two possible cases: Since we can always divide the equation Du = λ k u by a nonzero constant, we can take a = −1 in (i) and similarly we can normalize so that b = −1 in (ii). These two cases are discussed in Sections 5 and 6 respectively. In the next section we treat in details the case d = 1.

One dimensional orthogonal polynomials
The one dimensional case is the simplest case with no compatibility conditions. Our goal in this section is to provide a detailed study in this simplest situation and to show how the procedure works (i.e. how we determine the weight function W (x) and the set V ). The discrete orthogonal polynomials of one variable will also serve as building blocks in the higher-dimensional cases.
Historically, orthogonal polynomials satisfying difference equations of one variable were studied early in last century, we refer to [1,6] for references. It is well known that there are four families of such orthogonal polynomials of one variable, namely, Hahn polynomials, Meixner polynomials, Krawtchouk polynomials, and Charlier polynomials. However, to our great surprise, another family of orthogonal polynomials of discrete variable shows up in our study.
For simplicity of notation, we will drop all indices within this section, i.e. x = (x 1 ) ∈ R, E = E 1 , △ = △ 1 etc. The operator D will be of the form where (4.2) As we explained at the end of the previous section, we have two possible cases: (i) quadratic eigenvalue λ k , we can fix a = −1 (ii) linear eigenvalue λ k , we can fix a = 0, b = −1.
Hahn polynomials. Let (4.4) vanish for x = N + 1 but be positive for x ≤ N . Then one of α 1 and α 2 must be equal to −N − 1, say α 2 = −N − 1 and we can rewrite (4.4) as where we denoted β 1 = −α 3 − N − 1. Thus, up to a constant factor, the weight becomes The corresponding orthogonal polynomials on V 1 N = {0, 1, . . . , N } are the Hahn polynomials Q n (x; α 1 , β 1 , N ) given by In order to have a positive weight at all x ∈ V 1 N the right-hand side of (4.5) must be positive for x = 1, 2, . . . , N . Putting x = 1 and x = N we see that (α 1 +1)(β 1 +N ) > 0 and (α 1 + N )(β 1 + 1) > 0. Starting with β 1 > −1 or β 1 < −1 shows that there are two solutions, Conversely, it is clear that if α 1 and β 1 satisfy (i) or (ii) the weight function W (x) will have a constant sign on V 1 N . Hahn polynomials satisfy the following orthogonal relation The difference operator D takes the form Hahn-type polynomials on N 0 . Assume now that the right hand side of (4.4) does not vanish for any x ∈ N. Then it must be positive for all x ∈ N. In order to have a positive function W (x) and to have convergent series x∈N0 x n W (x) for some n ∈ N 0 , the parameter α 3 must be positive and α 1 , α 2 must satisfy one of the following (i) α 1 and α 2 are both positive, or there exists a negative integer κ such that α 1 , α 2 ∈ (κ, κ + 1); (ii) α 2 =ᾱ 1 ∈ C \ R. For such α 1 , α 2 , α 3 we can write the weight W (x) as The series x≥0 x n W (x) converges absolutely for n < α 3 − α 1 − α 2 − 1, i.e. the corresponding polynomials will be orthogonal on N 0 up to a given degree. It seems that these polynomials have not appeared in the literature before. Therefore, we derive below their basic properties: explicit formula in terms of hypergeometric functions, the three term recurrence formula, etc.
Proof. We look for a solution of (4.11) of the form Notice that (4.14) Using these formulas we can easily calculate D(R n ): where in the last formula we assumed that a n+1 = 0. Thus the difference equation (4.11) leads to the recursive relation whose solution with a 0 = 1 is This shows that the polynomials R n (x) defined by (4.10) satisfy (4.11). The orthogonality follows from Remark 3.4.
From the explicit formula we can derive other properties of these orthogonal polynomials.
Proof. Using again (4.14) we see that where a k are given by (4.15). It is now a simple matter to compare it with a n R n+1 (x) − (a n + c n )R n (x) + c n R n−1 (x) and verify that the coefficients of m k (x) agree for k = 0, 1, . . . , n + 1.

(4.18)
Proof. Using the three-term relation, a standard computation shows that R n , R n = c n c n−1 . . . c 1 a n−1 a n−2 . . . a 0 R 0 , R 0 , The quantity R 0 , R 0 = 1, 1 is a 2 F 1 evaluated at 1, which can be computed in closed form using Gauss' formula The statement follows from the above relations, (4.17) and the formula Remark 4.4. Notice that although R n , R m converges only for finitely many (n, m), formula (4.10) defines polynomials for every n ∈ N 0 . The operator D given in (4.12) is admissible according to Definition 3.3. Moreover, equations (4.11) and (4.16) hold for every n ∈ N 0 . In particular, this implies that the polynomials R n (x) provide a solution to the discrete-discrete version of the bispectral problem discussed in [3].  Assume now that l = 1. Then, we can rewrite (4.19) as where we have put c = (l − 1)/l = 0 and d = (1 + s − l)/(l − 1). We now have two possibilities depending on the sign of c.
Krawtchouk polynomials. If c < 0 then (4.21) can be rewritten as from which it follows that (−d) must be a positive integer (otherwise, W (x) will become negative at some point). If we put d = −N − 1 and p = c/(c − 1) ∈ (0, 1) we obtain the Krawtchouk polynomials The difference operator takes the form Meixner polynomials. If c > 0 then W (1) > 0 implies that β := d + 1 > 0. In order to have convergent series, we need c < 1. This leads to the Meixner polynomials depending on two parameters 0 < c < 1 and β > 0, orthogonal on N 0 with respect to the weight W (x) = c x (β) x /x!. The difference operator is The next proposition gives the general solution to the compatibility conditions (2.12) Proposition 5.1. For a = −1 there are 2 solutions to the compatibility conditions. Solution 1: where b, r, and {l i } d i=1 are free parameters.
Solution 2: Proof. Let us take i = j and let k = i, j. Comparing the coefficients of x 3 i x j x k and x 3 i x j on both sides of (2.12a) we get: − (l j,k + l j,i + 2) = −2(l j,i + 1), (5.4) and (5.5) equation (5.4) simply means that for i = j, l j,i is independent of i, i.e. we can put l j := l j,i for j = i. Using this in (5.5) we obtain that Our formulas simplify as follows With the above formulas it is very easy to solve completely (2.12a). Indeed, after canceling the common factor Let us denote r i = l i,i − l i . Then the first factor in the last formula is simply r i − r j . From (5.8) it follows that we have 2 possibilities: r i = r j for all i = j, i.e. r := r i is independent of i, or b = d k=1 l k −1. The first leads to Solution 1, and the second leads to Solution 2. It is immediate to see that the functions defined by (5.2) or (5.3) satisfy (2.12b).
Next, we want to determine the possible weights and orthogonal polynomials. It is easy to see that there is no weight corresponding to Solution 2. Otherwise, from (5.3) and (2.9) we will have Then, for x = e i it follows that l i < 0, i.e we can put −l i = r i > 0 and therefore, up to a constant factor, W (x) is given by but then even 1, 1 will diverge. Let us now concentrate on Solution 1, given by (5.2). Equations (5.2) and (2.9) give . (5.14) We have several possible sub-cases, depending on whether the right hand side in equation (5.14) can vanish or not. We divide them first in two big subclasses: 5.1. Case 1: b + r ∈ N 0 . Then |x| − b − r − 1 does not vanish. We have two possibilities: 5.1.1. For some i ∈ {1, 2 . . . , d}, (−σ i ) ∈ N. We show below that this implies that σ j > −1 for all j ∈ {1, 2 . . . , d}. Take j = i, arbitrary k ∈ N and n ∈ N 0 . From (5.14) we see that The second ratio is positive for large n, and therefore k + σ i > 0 for every k ∈ N, i.e. σ i > −1. Equation (5.15) for n = 0 implies that for every k ∈ N. But then for every j ∈ {1, 2, . . . , d} we must have showing that σ j > −1.
Since σ j > −1 for all j and b + r / ∈ N 0 we must have V = N d 0 . Up to a constant, the weight W (x) is given by Changing the parameters the weight takes the form and the free parameters are {σ i }, β and γ. The difference operator becomes and λ n = n(−n + γ − β − |σ| − d).
where α 2,j and α 3,j are given by satisfy the difference equation Dφ ν = λ |ν| φ ν and the orthogonal relation Proof. Notice that α 3,j − α 2,j is independent of x, and therefore from (4.10) one can immediately see that φ ν is indeed a polynomial of x of total degree |ν|. From Theorem 3.7 we know that for every k ∈ N 0 , the equation Du = λ k u has r k = k+d−1 k = dim(Π d k /Π d k−1 ) linearly independent solutions. Therefore, it is enough to prove that (5.18) holds.
In (5.18) we will first sum with respect to x d , then with respect to x d−1 , and so on. Writing (β + 1) |x| = (β + 1 + |X d−1 |) x d (β + 1) |X d−1 | and extracting only the terms depending on x d we get the sum . Using now (4.18) and the fact that α 3,d − α 2,d = γ − β is independent of x, we see that where B ν d is a constant (independent of x), whose value can be extracted from (4.18) .
Using the fact that 20) we see that . This shows, in particular, that where we have used the identity Equation (5.21) shows that the remaining d − 1 fold sums of R ν , R µ have exactly the same structure as that of the original d fold sums with β and γ replaced by β + ν d and γ − σ d − ν d − 1, respectively. In other words, it shows that we can use induction to complete the proof. For A ν we obtain the following formula: . (5.23) Using several times (5.20) one can rewrite the right-hand side as in (5.19).
From this we see that the parameters β and r must satisfy one of the following conditions: (i) β > −1 and r > 0; (ii) β < −|l| and r < −|l| + 1. Recall that b + r / ∈ N 0 , but everything will hold even if b + r = −β − 1 ∈ N 0 as long as β + |l| < 0 (i.e. if (ii) holds), because (β + 1) |x| will not vanish for x ∈ V d l . In the following we will use the notations y j and Y j defined in (5.17). For example, |L j | = l 1 + · · · + l j . If ν, l ∈ N d 0 , then ν ≤ l means ν i ≤ l i for 1 ≤ i ≤ d. Recall that the Hahn polynomials are denoted by Q n (x; α 1 , β 1 , N ).
(α 1,j + 1) νj Q νj (x j ; α 1,j , α 2,j , l j ), (5.25) where α 1,j and α 2,j are given by satisfy the difference equation Dφ ν = λ |ν| φ ν and the orthogonal relation where the normalization constant is given by By the definition of α 2,d and expanding the Pochhammer symbols gives Consequently, the weight function W d−1 becomes Apart from a constant multiple, the weight function W d−1 has exactly the same structure of W (x) with β replaced by β + ν d , r replaced by r + ν d , l replaced by L d−1 , respectively, and one variable less. Hence, we can use induction to complete the proof. For B ν we obtain the following formula gives (5.27).
The above ratio must be positive for all x ∈ V d N . In particular, for x = e i and x = N e i we see that From this it follows easily that the parameters {σ i } d+1 i=1 satisfy one of the following conditions: (i) σ i > −1 for i = 1, 2, . . . , d + 1; (ii) σ i < −N for i = 1, 2, . . . , d + 1.
The weight function takes the form with {σ i } as free parameters. If σ i > −1 for i = 1, 2, . . . , d + 1 or if σ i < −N for i = 1, 2, . . . , d + 1 but N is even then W (x) > 0 on V d N . If σ i < −N for i = 1, 2, . . . , d + 1 and if N is odd then W (x) < 0 on V d N , so one needs to change the sign in formula (5.29), in order to get a positive function.
The difference operator takes the form and the eigenvalues are λ n = −n(n + |σ| + d).
Remark 5.6. Notice that (5.25) in the parallelepiped case V d l gives a polynomial of total degree |ν|. We can use the generators G d l in Example 3.1 to write the same polynomial using only the monomials x µ with µ ∈ V d l . Similarly, in the case of Theorem 5.5 we can use the generators {(−x i ) li+1 : i ∈ S} to express Q ν in terms of the monomials x µ with µ ∈ V d N,S . 5.3. Summary. If a = −1 we have essentially one difference equation with coefficients given in (5.2). By specifying the free parameters, however, we end up with four different types of solutions given in Theorems 5.2 -5.5, respectively.
As an example, let us consider the case d = 2. The weight functions and the corresponding sets on which they live are listed below: In the last case, one can also exchange x 1 and x 2 to get another case. It should be mentioned that the equation (5.2) was considered in [12], but only the case (iii) was identified there.
First, notice that we cannot have γ i = 0. Indeed, if we assume that γ i = 0 then in order to have a self-adjoint operator we need β i = 0, which implies that l i,k = 0 for k = i. But then α i,k = 0 for k = i and therefore α k,i = 0 for k = i, which simply means that the operator D is independent of x i and E ±1 i . Thus we have an operator acting in a (d − 1) dimensional space, trivially extended, by adding x i as explained in Remark 3.6.
Below we assume that γ i = 0 for all i. Notice that γ i depends only on x i , so in (2.12b) we can cancel γ i γ j and we get and plugging these in (6.2) we get for i = j. Comparing the coefficients of x k in the last formula we see that These two equations have essentially two different types of solutions, which we discuss in two subsections.
Then the most general solution of the compatibility conditions (2.12) is given by x k + s (6.6c) where s and {l i } d i=1 are free parameters. Proof. From equation (6.5) it follows that l j,i = l j,j − 1. Denote l j := l j,j − 1. Then we have l j,i = l j for all i = j. Using now (6.3) we see that l j s i = l i s j . Thus s i = l i s, which leads to formulas (6.6).
Conversely, it is straightforward to see that if we define α i,j , β i and γ i as in (6.6a), (6.6c) and (6.6d), then the compatibility conditions in Corollary 2.3 are satisfied, i.e. the above formulas give the most general solution in the case l i,j = 0.
Below we determine the weight functions and the corresponding orthogonal polynomials. For every i we have which shows that l i must have the same signs. There are two possible cases.
6.1.1. l i = −p i < 0 for all i. In this case, If the denominator does not vanish, then for |x| large the second ratio will be negative and therefore we must have 1 − |p| < 0, i.e. |p| > 1. But then if we denote c i = p i /(|p| − 1) we will have |c| > 1 and up to a constant factor the weight is which leads to divergent series. Hence the only possibility here is (−s) = N ∈ N 0 . This forces |p| < 1 and The difference operator is then and the eigenvalues are λ n = −n. The orthogonal polynomials are the Krawtchouk polynomials on V d N . Recall that Krawtchouk polynomial in one variable is denoted by K n (x; p, N ).
The orthogonality of K ν (x; p, N ) follows from (5.31) under the limit. We refer the reader to [5] for other properties and applications of Krawtchouk polynomials. There is an analogous relation between the Meixner polynomials (6.11) of several variables and the orthogonal polynomials on the parallelepiped in (5.25), at least when parameters c have rational values. Indeed, let φ ν (x; β, r, l) denote the polynomials defined in (5.25), which are orthogonal with respect to the weight function in (5.24). In these polynomials we set β = s − 1 and Using induction, it is easy to see that l j /c j = l j−1 /c j−1 , which then implies that l j = (|l| + r)c j , or l j = c j 1 − |c| r, 1 ≤ j ≤ d.
If c j are rational numbers, we can choose r = (1 − |c|)N , N ∈ N 0 , for certain N , so that l j are integers. Upon taking r → ∞, which shows l j → ∞, it follows from (4.26) that lim l→∞ φ ν (x; s − 1, r, l) = M ν (x; s, c). (6.13) 6.2. Case 2: l i,j = 0 for some 1 ≤ i = j ≤ d. Recall that α i,j = 0 for i = j if and only if α j,i = 0. This essentially means that l i,j and l j,i are simultaneously zero or nonzero for i = j.
Lemma 6.5. If j = k and l j,k = 0 then for every i = j, k we have l i,j = l j,i = 0 or l i,k = l k,i = 0.
Proof. Follows immediately from (6.4). Proof. First we show that for m ∈ I and k ∈ J we have l m,k = l k,m = 0. Indeed if m = i, then l k,i = 0 by the definition of J. If m = i then l m,i = 0, but l k,i = 0. The previous lemma shows that l k,m = l m,k = 0. Thus if m ∈ I and k ∈ J we see at once that: • α m,k = α k,m = 0 • β m contains only the variables from I • β k contains only the variables from J. The decomposition of D follows immediately from the above observations. Remark 6.7. Theorem 6.6 says that if l i,j = 0 for some i = j, then the operator D splits as a sum of 2 operators of independent variables. Conversely let D I be an admissible operator in the variables x ′ = {x m : m ∈ I} and D J be an admissible operator in x ′′ = {x k : k ∈ J} with I ∩ J = ∅ and I ∪ J = {1, 2, . . . , d}. Let us denote by p I µ (x ′ ) and p J ν (x ′′ ) the polynomials satisfying D I (p I µ ) = −|µ|p I µ and D J (p J ν ) = −|ν|p J ν . Then if put D = D I + D J and p µ,ν (x) = p I µ (x ′ )p J ν (x ′′ ), where x = (x ′ , x ′′ ) ∈ R d , we have Dp µ,ν (x) = −(|µ| + |ν|)p µ,ν (x).
Thus, in this case, the eigenfunctions of the difference operator are product of orthogonal polynomials of fewer variables, which satisfy difference equations of lower dimension with linear eigenvalues. For d = 1 these are the polynomials of Charlier, Krawtchouk and Meixner given in Section 4.2. In higher dimensions, there are also the Krawtchouk polynomials of several variables given in Theorem 6.2 and the Meixner polynomials of several variables given in Theorem 6.3.
Clearly there are many product polynomials of this type and the number increases drastically as the dimension grows. As an example, we list all cases for d = 3 below. To list the different types we use the abbreviation of C, K, M for Charlier, Meixner, and Krawtchouk polynomial of one variable, respectively, and use K 2 and M 2 to denote the Krawtchouk polynomials of two variables on V 2 N and Meixner polynomials of two variables on N 2 0 . A product polynomial is denoted by its components. As an example, CCM stands for a product of the type Charlier-Charlier-Meixner.