Coupling coefficients of $su_q(1,1)$ and multivariate $q$-Racah polynomials

Gasper&Rahman's multivariate $q$-Racah polynomials are shown to arise as connection coefficients between families of multivariate $q$-Hahn or $q$-Jacobi polynomials. The families of $q$-Hahn polynomials are constructed as nested Clebsch--Gordan coefficients for the positive-discrete series representations of the quantum algebra $su_q(1,1)$. This gives an interpretation of the multivariate $q$-Racah polynomials in terms of $3nj$ symbols. It is shown that the families of $q$-Hahn polynomials also arise in wavefunctions of $q$-deformed quantum Calogero--Gaudin superintegrable systems.


Introduction
This paper shows that Gasper & Rahman's multivariate q-Racah polynomials arise as the connection coefficients between two families of multivariate q-Hahn or q-Jacobi polynomials. The two families of q-Hahn polynomials are constructed as nested Clebsch-Gordan coefficients for the positive-discrete series representations of the quantum algebra su q (1,1). This result gives an algebraic interpretation of the multivariate q-Racah polynomials as recoupling coefficients, or 3n j-symbols, of su q (1,1). It is also shown that the families of q-Hahn polynomials arise in wavefunctions of q-deformed quantum Calogero-Gaudin superintegrable systems of arbitrary dimension.
The multivariate q-Racah polynomials considered in this paper were originally introduced by Gasper and Rahman in [7] as q-analogs of the multivariate Racah polynomials defined by Tratnik in [34,35]. These q-Racah polynomials sit at the top of a hierarchy of orthogonal polynomials that extends the Askey scheme of (univariate) qorthogonal polynomials; Tratnik's Racah polynomials and their descendants similarly generalize the Askey scheme at q = 1. These two hierarchies will be referred to as the Gasper-Rahman and Tratnik schemes, respectively. The Gasper-Rahman scheme of multivariate q-orthogonal polynomials should be distinguished from the other multivariate extension of the Askey scheme based on root systems, which includes the Macdonald-Koornwinder polynomials [25] and the q-Racah polynomials defined by van Diejen and Stokman [36].
Like the families of univariate polynomials from the Askey scheme, the polynomials of the Gasper-Rahman and the Tratnik schemes are bispectral. Indeed, as shown by Iliev in [17], and by Geronimo and Iliev in [14], these polynomials simultaneously diagonalize a pair of commutative algebras of operators that act on the degrees and on the variables of the polynomials, respectively. The bispectral property is a key element in the link between these families of polynomials, superintegrable systems, recoupling of algebra representations, and connection coefficients of multivariate orthogonal polynomials. Recall that a quantum system with d degrees of freedom governed by a Hamiltonian H is deemed maximally superintegrable if it admits 2d − 1 algebraically independent symmetry operators, including H itself, that commute with the Hamiltonian [26].
For the univariate Racah polynomials, one has the following picture [11]. First, upon considering the 3-fold tensor product representations of su (1,1), one finds that the two intermediate Casimir operators associated to adjacent pairs of representations in the tensor product satisfy the (rank one) Racah algebra, which is also the algebra generated by the two operators involved in the bispectral property of the univariate Racah polynomials. This leads to the identification of the Racah polynomials as 6 j (or Racah) coefficients of su (1,1), which are the transition coefficients between the two eigenbases corresponding to the diagonalization of the intermediate Casimir operators. Second, if one chooses the three representations being tensored to belong to the positivediscrete series, the total Casimir operator for the 3-fold tensor product representation can be identified with the Hamiltonian of the so-called generic superintegrable system on the 2-sphere, and the intermediate Casimir operators correspond to its symmetries. Finally, one obtains the interpretation of the univariate Racah polynomials as connection coefficients between two families of 2-variable Jacobi polynomials that arise as wavefunctions of the superintegrable Hamiltonian. For a review of the connection between the Askey scheme and superintegrable systems, see [21]. For a review of the approach just described, the reader can also consult [10,12].
The picture described above involving the one-variable Racah polynomials has recently been fully generalized to Tratnik's multivariate Racah polynomials. In [20], Iliev and Xu have shown that these polynomials arise as connection coefficients between 1. Basics of su q (1,1) This section provides the necessary background material on the quantum algebra su q (1,1). In particular, the coproduct and the intermediate Casimir operators are introduced, and the representations of the positive-discrete series are defined.

su q (1, 1), tensor products and Casimir operators
Let q be a real number such that 0 < q < 1. The quantum algebra su q (1, 1) has three generators A 0 , A ± that satisfy the defining relations recovers the defining relations of su q (1, 1) in their usual presentation, that is The Casimir operator Γ, which commutes with all generators, has the expression The coproduct map ∆ : su q (1, 1) → su q (1, 1) ⊗ su q (1, 1) is defined as The coproduct ∆ can be iterated to obtain embeddings of su q (1, 1) into higher tensor powers. For a positive integer d, let ∆ (d) : su q (1, 1) → su q (1, 1) ⊗d be defined by where Id stands for the identity; one has ∆ (2) = ∆.
To each set [i; j], one can associate a realization of su q (1, 1) within su q (1, 1) ⊗d . Denoting by A (k) 0 , A (k) ± the generators of the k th factor of su q (1, 1) in su q (1, 1) ⊗d , the realization associated to the set S = [i; j] has for generators To each set S = [i; j], one can thus associate an intermediate Casimir operator Γ S defined as For a given value of d, the Casimir operator Γ [1;d] will be referred to as the full Casimir operator.

Representations of the positive-discrete series
Let α > 0 be a positive real number and let V (α) be an infinite-dimensional vector space with orthonormal basis e (α) n , where n is a non-negative integer. The space V (α) supports an irreducible representation of su q (1, 1) defined by the actions with 〈e i , e j 〉 = δ i j and where σ n is given by Note that σ (α) 0 = 0 and that when 0 < q < 1, one has σ (α) n > 0 for n ≥ 1. It follows that The su q (1, 1)-modules V (α) belong to the positive-discrete series. On V (α) , the Casimir operator (2) acts as a multiple of the identity. Indeed, one easily verifies using (2) and (6) that With the help of the nested coproduct defined in (3), one can define tensor product representations of su q (1, 1). Let α = (α 1 , α 2 , . . ., α d ) with α i > 0 be a d-dimensional multi-index and let W (α) be the d-fold tensor product As per (4), the space W (α) supports a representation of su q (1, 1) realized with the generators A . The space has a basis e (α) y defined by e (α) y = e (α 1 ) where y = (y 1 , . . ., y d ) is a multi-index of non-negative integers. The action of the generators A on the basis vectors (10) is easily obtained by combining (4) and (6). As a representation space, W (α) is reducible and has the following decomposition in irreducible components: When d = 2, a detailed proof of the decomposition (11) can be found in [33], see where γ(α) is given by (8); these eigenvalues have multiplicity m k .

Multivariate q-Hahn bases and the Clebsch-Gordan problem
In this section, two orthogonal bases of multivariate q-Hahn polynomials are constructed in the framework of the generalized Clebsch-Gordan problem for su q (1,1). The connection between the standard Clebsch-Gordan problem involving two-fold tensor product representations and the one-variable q-Hahn polynomials is first reviewed, and then generalized to the multifold tensor product case. Two related bases of multivariate q-Jacobi polynomials are also introduced through a limit.
Remark 2. Let us note that the functions C (α 1 ,α 2 ) n 1 ,n 2 (y 1 , y 2 ) are bispectral. Indeed, if we consider y and n such that y 1 + y 2 = n 1 + n 2 = N, then in view of (13), we can think of C (α 1 ,α 2 ) n 1 ,n 2 (y 1 , y 2 ) as an orthonormal q-Hahn polynomial of the variable y 1 , with index n 1 , multiplied by the square root of the weight. In addition to the spectral equation in y 1 stemming from the realization (18), they obey also a recurrence relation in the index n 1 . This property, which is well known (see [23,Section 14.6]), can be derived explicitly in the present context by considering the matrix element n 1 ,n 2 (y 1 , y 2 ), and by computing the action of q A (1) 0 on f (α 1 ,α 2 ) n 1 ,n 2 .

Nested Clebsch-Gordan coefficients and q-Hahn bases
We now consider the general d-fold tensor product representation W (α) defined in (9), and its decomposition in irreducible components (11). Since there are multiplicities in the decomposition, the eigenvalue problem for the total Casimir operator Γ [1;d] is degenerate. In the following, we will construct two bases associated to the diagonalization of two different sequences of intermediate Casimir operators. , the resulting basis will be orthogonal. Taking n = (n 1 , . . ., n d ), we define the orthonormal basis g (α) n of W (α) by the eigenvalue equations

First basis
where N k and A k are defined as We shall consider the expansion coefficients of the coupled basis g (α) n in the direct product basis e (α) y . Upon iterating the Clebsch-Gordan decomposition (12), one obtains the following result.
n (y) be defined by the expression where Y k = k i=1 y i , and where h n (x, α, β, N; q) is given by (14). These functions satisfy the orthogonality relations and arise in the expansion formulas Proof. One starts from the direct product basis e (α) y of W (α) . One can diagonalize the intermediate Casimir operator Γ [1;2] using the expansion (12). This leads to Upon using Remark 1 to identify f ;2] and Γ [1;3] simultaneously. One then obtains Using Remark 1 again to identify f [1;4] simultaneously. This leads to e (α) y = n 1 ,n 1 n 2 ,n 2 n 3 ,n 3 C (α 1 ,α 2 ) n 1 ,n 1 (y 1 , y 2 )C (2n 1 +α 1 +α 2 +1,α 3 ) n 2 ,n 2 (n 1 , y 3 ) and so on. Then one can use (13) to deduce that the terms in the sum above vanish unless Upon substituting the above condition in the expansion of e (α) y we obtain the second formula in (22) where the coefficients Ψ (α) n (y) are given in (20). Since both bases {e (α) y } and {g (α) n } are orthonormal, (21) and the first formula in (22) follow from the fact that the matrix (Ψ (α) n (y)) n,y is orthogonal.
Using (18), we can represent the action of Γ [1;k] as a difference operator in the variables y. Then, the equation combined with (19) shows that Ψ (α) n (y) are eigenfunctions of the difference operators Γ [1;k] with eigenvalues γ(2N k−1 + A k +k−1). The results in [17] imply that the functions Ψ (α) n (y) are also eigenfunctions of commuting operators acting on the variables n, and therefore are bispectral. The natural extension of the arguments in Remark 2 provides a Lie-interpretation of the corresponding bispectral algebras of partial difference operators.

A family of multivariate q-Hahn polynomials
The orthonormal functions Ψ (α) n (y) can be expressed in terms of the multivariate q-Hahn polynomials introduced by Gasper and Rahman in [7]. Indeed, one can write with ρ (α) (y) given by The normalization factor Λ (α) n has the expression and H (α) n (y) are the Gasper-Rahman multivariate q-Hahn polynomials defined as where h n (x, α, β, N; q) are the q-Hahn polynomials given in (15). These multivariate q-Hahn polynomials are a direct q-deformation of Karlin & McGregor's multivariate Hahn polynomials [22]. Upon fixing M ∈ N and taking n and n ′ such that N d = N ′ d = M, these polynomials satisfy the orthogonality relation Remark 3. Note that when dealing with the polynomials (25) [2;3] , Γ [2;4] , . . ., Γ [2;d] , Γ [1;d] }. These operators are self-adjoint on W (α) , and they commute with one another. Consequently, one can construct an orthonormal basis that simultaneously diagonalizes them. Taking m = (m 1 , . . ., m d ), the orthonormal basis u (α) m of W (α) is defined by the eigenvalue equations

Second basis Consider the sequence of operators {Γ
where Once again we consider the expansion coefficients of the coupled basis u (α) m in the direct product basis e (α) y . One has the following result.
where Y k = k i=2 y i and where h n (x, α, β, N; q) is given by (14). The function Ξ (α) m (y) satisfy the orthogonality relations and arise in the expansion formulas Proof. The proof is along the same lines as that of Proposition 1. Starting from the direct product basis e (α) y , one diagonalizes the operator Γ [2;3] by using the expansion (12). This leads to One can use Remark 1 to identify f Using the explicit expression (13) then yields (28).
The coefficients Ξ (α) m (y) can also be viewed as nested Clebsch-Gordan coefficients and satisfy bispectral equations.

Another family of multivariate q-Hahn polynomials
The orthogonal functions Ξ (α) m (y) can also be written in terms of a family multivariate orthogonal polynomials of q-Hahn type. One has indeed where ρ (α) (y) is given by (24). The normalization factor Ω (α) m is of the form and the polynomials G (α) m (y) have the expression where h n (x, α, β, N; q) are the q-Hahn polynomials (15). These polynomials are orthogonal with respect to the same measure as the Gasper-Rahman q-Hahn polynomials H (α) n (y) given by (25). Upon fixing L ∈ N and taking m and m ′ such that M d = M ′ d = L, the orthogonality relation for the polynomials G (α) m (y) reads Let us quickly recap the results obtained so far. We have used nested Clebsch-Gordan coefficients for multifold tensor product representations of su q (1, 1) to construct two bases of multivariate orthogonal functions Ψ (α) n (y) and Ξ (α) m (y) that diagonalize two commutative subalgebras of intermediate Casimir operators. From Ψ (α) n (y) and Ξ (α) m (y), we then constructed two families of multivariate q-Hahn polynomials H (α) n (y) and G (α) m (y) that are orthogonal with respect to the same measure ρ (α) (y) given in (24).

q-Jacobi bases
The orthogonal functions Ψ (α) n (y) and Ξ (α) m (y) defined in (20) and (28) are both non-zero if the condition Y d = N d = M d is satisfied. The corresponding families of multivariate q-Hahn polynomials H (α) n (y) and G (α) m (y) consequently satisfy the finite orthogonality relations (26) and (33). It is therefore meaningful to consider limits of these families of functions and polynomials as Y d = N d = M d goes to infinity. We shall consider the limits of the functions Ψ (α) n (y) and Ξ (α) m (y) and extract orthogonal polynomials from each of those limits. These results will prove useful later.
2.3.1. A first basis of q-Jacobi functions Let us first consider the limit of the functions Ψ (α) n (y) as Y d = N d goes to infinity. We find the following result.
The functions J (α) s (x) obey the orthogonality relation where the multi-index x runs over multi-indices of non-negative integers.
Proof. The orthogonality relation (37) follows from (35). The limit can be taken directly. The calculation is long, but otherwise straightforward. The following result relating the q-Hahn to the q-Jacobi polynomials is useful [6]:

Interbasis expansion coefficients and q-Racah polynomials
In this section, we shall consider the expansion coefficients between the bases Ψ (α) n (y) and Ξ (α) m (y) and show that they are expressed in terms of Gasper & Rahman's multivariate q-Racah polynomials. Given the interpretation of the functions Ψ (α) n (y) and Ξ (α) m (y) as nested Clebsch-Gordan coefficients of the d-fold tensor product representation W (α) , this will provide an interpretation of the multivariate q-Racah polynomials in terms of 3n j symbols for su q (1, 1).
Here also, the coefficients are independent of the coordinates x.

Connection coefficients between multivariate polynomials
The functions R (α) m (n) are connection coefficients between families of multivariate orthogonal polynomials. Indeed, it follows from (23), (31) and (45) that 1 where the sum is over the multi-indices n such that N d = M d and N d−1 = M d−1 . One can also use the orthogonality relation for H (α) n (y) to get the formula where the sum is restricted to all y such that Y d = N d = M d . Furthermore, from (39), (43) and (48) one finds that which is equivalent to As can be seen, the coefficients R (α) m (n) serve as connection coefficients between bases of multivariate q-Hahn or q-Jacobi polynomials. Since these bases are themselves orthogonal, it follows from elementary linear algebra that the functions R (α) m (n) satisfy the orthogonality relations These relations are meaningful when the indices n, n ′ , m and m ′ are such that

The d = 3 case: one-variable q-Racah polynomials
Let us now consider the expansion coefficients R (α) m (n) for d = 3. As we noted above, the coefficients R (α) m (n) vanish unless M 2 = N 2 and m 3 = n 3 , and are independent of the quantum numbers m 3 and n 3 . Throughout this subsection, we assume that m and n satisfy these conditions, and to simplify the notation, we will omit m 3 and n 3 when we display the coefficients, i.e. we will write simply R (α 1 ,α 2 ,α 3 ) m 1 ,m 2 (n 1 , n 2 ) instead of R (α 1 ,α 2 ,α 3 ) m 1 ,m 2 ,m 3 (n 1 , n 2 , n 3 ).
The calculation of the normalization factors is straightforward.
Let us note that our derivation of Proposition 5 is much simpler than the one presented in [4], which used different methods. Moreover, our construction has the advantage of unifying two of the main interpretations of the one-variable q-Racah polynomials: 1) their interpretation as connection coefficients for 2-variable q-Hahn or q-Jacobi polynomials, and 2) their interpretation as 6 j coefficients for positivediscrete series representations of su q (1, 1). In the next section, we shall give a new interpretation in connection with q-deformed quantum superintegrable systems.
Remark 7. We can use special values of the y variables in equation (49) to obtain other generating functions for the q-Racah polynomials. First, note that by using the qanalog of the Chu-Vandermonde identity, the q-Hahn polynomials in (15) reduce to simple products of q-Pochhammer symbols when x = −(α + 1) or x = N: h n (N, α, β, N; q) = (q −N ; q) n (q −n−β ; q) n q n(n+α+β+1) .

The general result
Proof. By induction; following the steps outlined in the previous subsection. The basic case d = 3 is proven in Proposition 5. Suppose that the result holds at level d − 1.
Consider the basis functions Ξ (α) m (y) at level d. One can write Then, one uses the induction hypothesis to develop the functions Ξ We have thus obtained the explicit expression (59) for the expansion coefficients R (α) m (n) between the q-Hahn bases Ψ (α) n (y) and Ξ (α) m (y) defined in (20) and (28). These are also the expansion coefficients between the q-Jacobi bases J (α) n (y) and Q (α) m (y) defined in (36) and (42). Moreover, since these coefficients are the overlaps between basis vectors corresponding to irreducible decompositions of the multifold tensor product representation W (α) of su q (1, 1), the coefficients R (α) m (n) can also be considered as particular 3n j-coefficients. Remark 8. Using (59), we show in the next subsection that the coefficients R (α) m (n) can be expressed in terms of multivariate q-Racah polynomials defined by Gasper and Rahman in [7]. Similarly to the d = 3 case discussed in Remark 7, we can use special values of the y variables in equation (49) to obtain different identities for these polynomials. However, if we fix the values of y 1 , . . ., y d , so that the q-Hahn polynomials reduce to products of q-Pochhammer symbols, there will be just one free variable left, which is not sufficient to characterize the q-Racah polynomials in the multivariate setting.

Multivariate q-Racah polynomials
The coefficients R (α) m (n) can be expressed in terms of the multivariate q-Racah polynomials introduced by Gasper and Rahman in [7]. In the s-variable case, these polynomials can be written as where y 0 = 0, y s+1 = M, and M ∈ N and β = (β 0 , β 1 , . . ., β s+1 ) are parameters. It is not hard to see that Z ℓ (y; β, M; q) is a polynomial of total degree L s in the variables z k = q −y k + β k q y k . Moreover, they are orthogonal on the simplex with respect to the weight function The parametrization in [7] can be obtained by taking a 1 = q β 1 , a k = q β k −β k−1 for k = 2, . . ., s + 1 and b = q β 1 −β 0 −1 . The square of the norm is Upon comparing (60) with (59), it is seen that the coefficients (59) will be proportional to the polynomials (60) if one takes s = d − 2 and The expansion coefficients (59) can be expressed in terms of the multivariate q-Racah polynomials (60) as follows: where the connection between the original variables and parameters n, m, α and y, ℓ, β, M is provided by (63).
Using formulas (61) Since the right-hand side of (70) is invariant under the involution (66), the last identity combined with Proposition 7 shows that q-Racah polynomials, viewed as polynomials in the indices, are also orthogonal with respect to an appropriate multivariate q-Racah weight (61). More precisely, the following statement holds.
The last equation combined with (64) and (65) relates the dual orthogonality relation (50b) to the orthogonality of the Gasper-Rahman q-Racah polynomials.

q-Deformed Calogero-Gaudin systems
In this section, it is shown that the bases of multivariate functions Ψ (α) n (y) and Ξ (α) m (y) constructed in Section 2 are in fact wavefunctions for quantum q-deformed Calogero-Gaudin superintegrable systems.
Let us return to the realization (18) of the quantum algebra su q (1, 1); i.e. we take q A (i) 0 = q y i +(α i +1)/2 , where T ± y i f (y i ) = f (y i ±1) is the discrete shift operator in the variable y i , and where σ (α) n is given by (7). Following the coproduct construction (4), one has su q (1, 1) realizations su S q (1, 1) associated to each set S = [i; j] acting on the variables y i , . . ., y j . These realizations read For a given value of d, the "full" su q (1, 1) realization corresponds to the set The Hamiltonian (73) corresponds to a q-deformed quantum Gaudin-Calogero system in (d − 1) dimension. These systems have been discussed in [28] in particular. Their integrability was shown, and the eigenvalues and a set of eigenvectors were obtained. In the present approach, it is clear that the Hamiltonian (73) is in fact superintegrable. Indeed, the elements of the two commutative subalgebras 〈Γ [1;2] (19) and (27). The multivariate q-Racah polynomials then correspond to the connection coefficients between two bases for the eigenstates of the quantum Calogero-Gaudin model (73).