Svante Janson

A divergent generating function that can be summed and analysed analytically

Keywords: recurrence, divergent generating function, Borel summation, Bessel functions, Lommel polynomials, holonomic function

We study a recurrence relation, originating in combinatorial problems, where the generating function, as a formal power series, satisfies a differential equation that can be solved in a suitable domain; this yields an analytic function in a domain, but the solution is singular at the origin and the generating function has radius of convergence 0. Nevertheless, the solution to the recurrence can be obtained from the analytic solution by finding an asymptotic series expansion. Conversely, the analytic solution can be obtained by summing the generating function by the Borel summation method. This is an explicit example, which we study detail, of a behaviour known to be typical for a large class of holonomic functions. We also express the solution using Bessel functions and Lommel polynomials.

Introduction

Consider the recurrence relation

a n+1 + a n-1 = (αn + β)a n , n ≥ 1, (1.1)
where α and β are given real (or possibly complex) numbers, and with initial values a 0 = 0, a 1 = 1, say.

We assume α = 0. (Otherwise, this is a linear recurrence with constant coefficients and the solution is well known. See e.g. Flajolet and Sedgewick (2008, Section IV.5.1).) Recurrence relations of this type appear in several contexts, see e.g. Renlund (2009+), Parviainen and Renlund (2010+) and the sequences A053983, A053984, A058797, A058798, A058799 in [START_REF] Sloane | The On-Line Encyclopedia of Integer Sequences[END_REF]; further examples from [START_REF] Sloane | The On-Line Encyclopedia of Integer Sequences[END_REF] are given in Table 5.

The recurrence (1.1) is easily solved using Bessel functions, see the references above; for completeness we give this solution in Appendix A. The present paper originates in an attempt to understand this solution using generating functions.

Define thus the generating function

A(z) := ∞ n=0 a n z n . (1.2)
It is easily seen that (1.1) with a 0 = 0 and a 1 = 1 is equivalent to the differential equation (cf. the general formulas in Section 3) z -1 A(z) + zA(z) = αzA (z) + βA(z) + 1 (1.3) or αz 2 A (z) = (1 + z 2 -βz)A(z) -z.

(1.4)

This differential equation is singular at 0. For z = 0 it can be written

A (z) -α -1 z -2 -βα -1 z -1 + α -1 A(z) + α -1 z -1 = 0.
(1.5)

We write, for convenience, ν := β/α and find the integrating factor H(z) := exp(α -1 z -1 + ν log zα -1 z) and the solution, in any given simply connected domain not containing 0, where we may choose any z 0 in the domain and C is an arbitrary constant,

A(z) = α -1 H(z) z0 z w -1 H(w) dw + C H(z)
= α -1 z -ν e -α -1 /z+α -1 z z0 z w ν-1 e α -1 /w-α -1 w dw + Cz -ν e -α -1 /z+α -1 z .

(1.6)

The standard procedure to find the coefficients from the generating function (1.2) is to find the Taylor coefficients at 0: a n = D n A(0)/n!. (See, e.g., [START_REF] Flajolet | Analytic Combinatorics[END_REF] where this is described and developed in detail.) But here this breaks down because A(z) is (in general) singular at 0. In fact, except in exceptional cases, the coefficients a n increase so rapidly that the sum (1.2) has radius of convergence 0; i.e., the sum (1.2) diverges for every complex z = 0, and there is no analytic function A(z) in any neighbourhood of 0 that satisfies (1.5). (See Remark A.4.) Hence, (1.2) should be interpreted as a formal power series, and in this sense (1.3) and (1.4) hold and determine A(z).

Nevertheless, we have in (1.6) found a solution (or, more precisely, a family of solutions) A(z) to (1.3) that is a function of a real or complex variable, and it is natural to ask in what way this solution represents the sequence (a n ) n , and if one can recover (a n ) n from A(z). We shall see that indeed there is a strong connection (at least if we choose the right solution, or restrict ourselves to a suitable domain): A(z) has an asymptotic expansion A(z) ∼ n a n z n as z → 0. (Recall that this means that A(z) = N n=0 a n z n + O(z N +1 ) as z → 0 for every N ≥ 0, but does not imply convergence of the infinite sum, see e.g., Flajolet and Sedgewick (2008, Section A.2).) Further, we will find a n explicitly from this expansion. Moreover, although A(z) is not analytic at 0, there is a solution A(z) of (1.3) that is infinitely differentiable on the entire real line, including 0, and the derivatives at 0 are a n n! as expected. We shall also see that although the sum n a n z n diverges, it can be summed by (a version of) the Borel summation method [START_REF] Hardy | Divergent Series[END_REF] for every z ∈ C \ (0, ∞). We can thus define the generating function A(z) in this domain by (1.2), interpreted as a (generalized) Borel sum, and this function A(z) satisfies (1.3)-(1.5), and thus (1.6) for some choice of z 0 and C. (In fact, we can take z 0 = 0 with some care in the path of integration in (1.6), and then C = 0.)

We shall prove these results for a somewhat larger class of recurrence relations in Sections 3 and 4. We then use this and (1.6) to find a simple formula for a n in Section 5. The appendices give another formula for a n using Bessel functions, and connect the two solutions.

As mentioned in Flajolet and Sedgewick (2008, Note A.10), the possibility to define a divergent formal power series as a function having the right asymptotic expansion goes back to [START_REF] Euler | De seriebus divergentibus[END_REF]; he considered a simpler example that we will treat in detail in Section 2 before studying the general case.

The generating function studied in this paper is an example of the class of holonomic functions, see Flajolet and Sedgewick (2008, Appendix B.4) for an introduction. Related and much more general results can be found in the books by [START_REF] Balser | From Divergent Power Series to Analytic Functions[END_REF][START_REF] Balser | Formal Power Series and Linear Systems of Meromorphic Ordinary Differential Equations[END_REF]; see further, e.g., [START_REF] Balser | Multisummability of formal power series solutions of linear ordinary differential equations[END_REF], [START_REF] Braaksma | Multisummability and Stokes multipliers of linear meromorphic differential equations[END_REF], [START_REF] Martinet | Elementary acceleration and multisummability[END_REF], [START_REF] Malgrange | Fonctions multisommables[END_REF][START_REF] Malgrange | Fonctions multisommables[END_REF]Ramis (1993).

Remark 1.1. By a theorem of [START_REF] Borel | Sur quelques points de la théorie des fonctions[END_REF], see also Carleman (1926, Ch. V), given any sequence (a n ), there exists a C ∞ function on R with these numbers as Taylor coefficients at 0, and thus the asymptotic expansion n a n z n there; moreover, we may choose the function so that it is analytic in, e.g., a given sector D in the complex plane, with the given asymptotic expansion as z → 0 in D. Hence, the existence of a function (and, indeed, infinitely many functions) representing a given sequence by an asymptotic expansion is well-known. The special feature in the problem discussed in this paper is the possibility to find such a function explicitly and to use it to find a formula for a n .

Acknowledgement. This research was partly done during a visit to the programme "Combinatorics and Statistical Physics" at the Erwin Schrödinger International Institute for Mathematical Physics in Vienna, 2008. The paper was presented at the colloquium for Philippe Flajolet's 60th birthday in Paris, December 2008. I am grateful both to several participants at these meetings, including Philippe Flajolet, and to Henrik Renlund for interesting comments and references.

A simple example

Let us first study a simpler recurrence where we easily can see in detail what happens:

b n = nb n-1 , n ≥ 1, with b 0 = 1. (2.1)
The solution is obviously b n = n!, but let us ignore this fact and try to find b n through the generating function B(z) = ∞ n=0 b n z n . Since the series diverges for any complex number z = 0, we regard B(z) as a formal power series. We have

∞ n=1 nb n-1 z n = ∞ n=0 (n + 1)b n z n+1 = ∞ n=0 nb n z n+1 + ∞ n=0 b n z n+1 = z 2 B (z) + zB(z), so (2.1) is equivalent to B(z) -1 = z 2 B (z) + zB(z) (2.2) or B (z) = (z -2 -z -1 )B(z) -z -2 . (2.3)
Now forget that the series defining B diverges, and let us regard (2.3) as a differential equation for a function B(z) of a real or complex variable. We find an integrating factor H(z) = e h(z) where h (z) = -z -2 + z -1 , i.e., we can choose h(z) = z -1 + log z and H(z) = ze 1/z ; this yields (for z = 0)

B(z)H(z) = B (z)H(z) + B(z)H (z) = -z -2 H(z) = -z -1 e 1/z ,
and thus, for any fixed z 0 = 0 and some constant C, B(z) = z -1 e -1/z z0 z w -1 e 1/w dw + Cz -1 e -1/z .

(2.4) (We can choose z 0 freely; a change in z 0 is compensated by a change in C.)

The formula (2.4) defines B(z) as an analytic function in any simply connected domain that does not contain 0, for example in the plane cut along the negative real axis or along any other ray from 0. Note that B cannot be defined as an analytic function in C \ {0}, since w -1 e 1/w dw = 2πi = 0 for a circle around 0 (by the power series expansion of e 1/w ).

Let us now see what happens as z → 0. Consider first the case when z < 0, i.e., z is real and negative. Note that e 1/z → 0 and e -1/z → ∞ rapidly as z 0. Hence, 0 z w -1 e 1/w dw converges, so we may choose z 0 = 0 in (2.4). Let us first also choose C = 0, and denote the resulting solution by B 0 (z); thus, using the change of variable t = 1/z -1/w, and thus w = z/(1 -zt),

B 0 (z) = z -1 e -1/z 0 z w -1 e 1/w dw = z -1 0 z e 1/w-1/z w -1 dw = ∞ 0 e -t 1 -zt dt, z < 0. (2.5)
By dominated convergence in the last integral, B 0 (z) → 0 as z 0. We made here a specific choice of the integration constant C in (2.4), but we now see that we really had no alternative; any other choice would give a solution B(z) such that |B(z)| → ∞ as z 0. Hence B 0 is the unique solution of (2.2), or (2.3), on (-∞, 0) that is continuous up to 0.

Moreover, using

(1 -zt) -1 = N 0 z n t n + O(z N +1 t N +1 ), (2.6)
in (2.5) and integrating, we find the asymptotic expansion

B 0 (z) ∼ ∞ n=0 n! z n = ∞ n=0 b n z n , z 0.
(2.7)

Alternatively, we can differentiate (2.5) repeatedly under the integral sign (using dominated convergence) and obtain, for any n ≥ 0,

D n B 0 (z) n! = ∞ 0 t n e -t (1 -zt) n+1 dt → ∞ 0 t n e -t dt = n!, as z 0.
(2.8)

Remark 2.1. Note that the asymptotic expansion (2.7) follows directly from (2.8) by Taylor's formula; however, an asymptotic expansion does not conversely imply convergence of the derivatives without further assumptions. (Consider e -1/x sin e 1/x 2 as x 0 as an example.) One situation where the converse holds is when the asymptotic expansion holds in a sector θ 1 < arg z < θ 2 in the complex plane; then Cauchy's estimates imply convergence of all derivatives in any strictly smaller sector θ 1 < arg z < θ 2 with θ 1 < θ 1 < θ 2 < θ 2 . (We use below this meaning of 'strictly smaller' or 'strictly contained in' for sectors.) More generally, or as a consequence, the asymptotic expansion can be formally differentiated in any strictly smaller sector. This shows one advantage of using complex variables, and will be used below.

Remark 2.2. So far, we have in this section essentially just repeated what Euler (1760) did in 1760 when he wanted to compute ∞ n=0 (-1) n n!. He took the sum to be, in our notation, B 0 (-1) = ∞ 0 e -t (1 + t) -1 dt, and then evaluated this integral numerically, by several methods, finding the approximate value 0.5963473621237 (where all but the last four digits are confirmed by Maple). See Hardy (1949, §2.4-2.6) for a modern presentation.

Next consider z > 0, so z 0 along the positive real axis. Then e -1/z decreases rapidly and is o(z N ) for every N , and the same holds for the term Cz -1 e -1/z in (2.4). We choose z 0 = 1 and again make the change of variable t = z -1 -w -1 , w = z/(1 -zt). This yields, for any C and N ,

B(z) = 1/z-1 0 e -t 1 -zt dt + Cz -1 e -1/z = 1/z-1 0 e -t 1 -zt dt + O(z N ), z > 0.
(2.9)

Assume 0 < z < 1/2. The integral from 1/(2z) to 1/z -1 is O(z N ). For 0 ≤ t ≤ 1/(2z), (2.6
) holds, and (2.9) yields the asymptotic expansion

B(z) ∼ ∞ n=0 n! z n = ∞ n=0 b n z n , z 0, (2.10)
just as (2.7) on the negative real axis, but holding for all solutions of (2.3) on the positive real axis. Again, we further can differentiate (2.9) under the integral sign (noting that we also get a term from the varying integration limit) and conclude that

D n B(z)/n! → n! = b n as z 0, for any solution B on (0, ∞).
Combining the solution B 0 (z) on (-∞, 0) and any solution B(z) on (0, ∞), and defining B(0) = b 0 = 1, we obtain a function B(z) defined on the entire real axis, which is infinitely differentiable except possibly at 0, and with all derivatives extending continuously also to 0; this implies (see the end of the proof of Theorem 3.2 for details) that B(z) is infinitely differentiable, and that it satisfies (2.2) on the entire real axis. We have thus found a family of solutions of (2.2) on R; we see that each such solution is infinitely differentiable on R and analytic except at 0.

These results for real z are easily extended to complex z. The arguments and results for z < 0 extend to any ray {z : arg z = θ} or sector {z : arg z ∈ (θ 1 , θ 2 )} strictly contained in the left half-plane, and the arguments and results for z > 0 extend to any ray {z : arg z = θ} or sector {z : arg z ∈ (θ 1 , θ 2 )} strictly contained in the right half-plane. (Note that we then can use the argument in Remark 2.1 as an alternative to show convergence of derivatives.) Moreover, the particular solution B 0 can be defined by (2.5) for all z ∈ C \ [0, ∞), and it is easily seen that (2.7) and (2.8) hold in this larger domain too. We omit the details, since we state and prove the more general Theorem 3.2 in the next section.

The general case

Consider, more generally, a linear recurrence

a n = K k=1 (α k n + β k )a n-k , n ≥ K, (3.1)
with given a 0 , . . . , a K-1 , where K ≥ 1 and α 1 , β 1 , . . . , α K , β K are given complex numbers with α 1 > 0.

Equivalently, the recurrence can be written as

a n = K k=1 α k (n -k) + βk a n-k , n ≥ K, (3.2)
with βk := β k + kα k . Clearly, (3.1) includes (1.1) and (2.1) as special cases.

Remark 3.1. The assumption α 1 > 0 is no essential restriction, and the results are easily extended to any α 1 = 0. In fact, we may as well assume α 1 = 1 when convenient; the general case α 1 = 0 is easily reduced to this by replacing

a n , α k , β k by α -n 1 a n , α k α -k 1 , β k α -k 1 .
Note, however, that unless α 1 > 0, this entails a rotation in the complex plane of the domains in Theorem 3.2 below; for example, if α 1 < 0, the intervals I + and I -have to be interchanged.

Let A(z) := ∞ n=0 a n z n , regarded as a formal power series. The recursive relation (3.2) is equivalent to the differential equation

A(z) = K k=1 α k z k+1 A (z) + K k=1 βk z k A(z) + r(z), (3.3)
where r(z) is a polynomial of degree less than K determined by the initial conditions. Introduce the polynomials p(z) := K k=1 α k z k-1 and q(z) := 1 -K k=1 βk z k ; note that p(0) = α 1 > 0 and q(0) = 1. Then (3.3) can also be written q(z)A(z) = z 2 p(z)A (z) + r(z) or

A (z) = q(z) z 2 p(z) A(z) - r(z) z 2 p(z) = R 1 (z)A(z) -R 2 (z), (3.4)
for two rational functions R 1 and R 2 . (iv) Any solution A(z) of (3.3) on I satisfies the asymptotic expansion A(z) ∼ ∞ n=0 a n z n as z → 0, i.e., for every N ,

A(z) = N n=0 a n z n + O(z N +1 ), z → 0. (3.5) (v) Every solution A(z) of (3.
3) or (3.4) on I + satisfies the asymptotic expansion (3.5) in I + , and

D n A(z)/n! → a n as z 0.
(vi) There is a unique solution A 0 (z) of (3.3) or (3.4) on I -that is bounded as z 0. This solution satisfies the asymptotic expansion (3.5), and D n A 0 (z)/n! → a n as z 0.

For complex z, these results extend as follows, for some δ > 0.

(If K = 1, or more generally α 2 = • • • = α K = 0, then δ can be chosen as ∞.) (vii)
The solution A 0 (z) extends to an analytic function in the domain {z : |z| < δ, z / ∈ [0, ∞)}; this solution to (3.3) and (3.4) satisfies the asymptotic expansion (3.5) and D n A 0 (z)/n! → a n as z → 0 in this domain.

(viii) For any θ 1 and θ 2 with -π/2 < θ 1 < θ 2 < π/2, there exist (non-unique) solutions A(z) of (3.3) and (3.4) in the domain {z : |z| < δ, arg(z) ∈ (θ 1 , θ 2 )}. Any such solution satisfies the asymptotic expansion (3.5) and D n A(z)/n! → a n as z → 0 in this domain.

(ix) For any θ 1 and θ 2 with π/2 < θ 1 < θ 2 < 3π/2, there exist solutions A(z) of (3.3) and (3.4) in the domain {z : |z| < δ, arg(z) ∈ (θ 1 , θ 2 )}. Exactly one of these solutions is bounded as z → 0, and this solution satisfies the asymptotic expansion (3.5) and D n A(z)/n! → a n as z → 0 in this domain.

Proof: Suppose for notational simplicity that α 1 = 1; the general case is easily reduced to this as explained in Remark 3.1. Then p(0) = q(0) = 1, and R 1 (z) = z -2 + γz -1 + g(z) for some constant γ and some rational function g that has poles only at the zeros of p, and thus not on I.

In a simply connected domain D in the complex plane where

z 2 p(z) = 0, a primitive function of -R 1 (z) is thus h(z) = z -1 -γ log z -G(z)
, where G = g. Consequently, the equation (3.4) has an integrating factor H(z) := e h(z) and the solutions in D are given by, for any fixed z 0 ∈ D and with a constant C ∈ C,

A(z) = H(z) -1 z0 z H(w)R 2 (w) dw + CH(z) -1 = z0 z e h(w)-h(z) R 2 (w) dw + Ce -h(z) = z0 z (z/w) γ e w -1 -z -1 -G(w)+G(z) R 2 (w) dw + Cz γ e -z -1 +G(z) .
(3.6)

Any choice of C gives an analytic solution in D, and all solutions in D are given by (3.6). The same is true if we consider solutions on an interval where z 2 p(z) = 0 (since any such interval is contained in a suitable domain D). Clearly, any solution of (3.3) or (3.4) on an interval or in a domain where z 2 p(z) = 0 is analytic. We write R 3 (z) = r(z)/p(z); thus R 2 (z) = z -2 R 3 (z) and R 3 is a rational function that is analytic at 0. We make in (3.6) the same change of variable t = 1/z -1/w as in Section 2; thus w = z/(1 -zt) and w -2 dw = dt, and we find

A(z) = 1/z-1/z0 0 (1 -zt) γ e -t-G(z/(1-zt))+G(z) R 3 z 1 -zt dt + Cz γ e -z -1 +G(z) .
(3.7) (We integrate here along any curve such that z/(1 -tz) ∈ D.) Consider first z ∈ I -. We choose D ⊃ I -and z 0 ∈ I -and integrate in (3.6) and (3.7) along parts of the real axis. As in Section 2, we can then let z 0 0, i.e., take z 0 = 0 in (3.6). It then follows from (3.6) that the solution is bounded as z 0 if, and only if, we choose C = 0. We denote this solution by A 0 (z) and find from (3.7)

A 0 (z) = ∞ 0 (1 -zt) γ e -t-G(z/(1-zt))+G(z) R 3 z 1 -zt dt. (3.8) Now consider z ∈ C \ [0, ∞). Then, for every t ≥ 0, |1 -zt| ≥ c(arg(z)), with c(θ) = | sin(θ)| if |θ| ≤ π/2, and c(θ) = 1 if |θ| ≥ π/2. Let ρ := min{|ζ| : p(ζ) = 0} (with ρ = ∞ if α 2 = • • • = α K = 0,
and thus p has no root). Then G(z), R 3 (z) and F (z) := R 3 (z)e -G(z) are analytic in {z : |z| < ρ}, and thus, for any ε > 0,

F (z/(1 -zt)) is analytic in D ε := {z : | arg(z)| > ε, |z| < c(ε)ρ} for every t ≥ 0.
If we take the factor e G(z) outside the integral in (3.8), we can thus differentiate the remaining integral any number of times under the integral sign, and the derivatives will, by induction, be linear combinations of integrals of the type

∞ 0 t (1 -zt) γ F (m) (z/(1 -zt
))e -t dt. By dominated convergence, each such integral converges as z → 0 in D ε , and it follows that every derivative D n A 0 (z) converges as z → 0 in D ε . Denote the limits by a * n , n ≥ 0. It follows by Taylor expansions that A 0 (z) has an asymptotic expansion (3.5) in D ε , with a n replaced by a * n , and that A 0 (z) has a similar asymptotic expansion, that can be obtained by formally differentiating the expansion for A 0 (z). Substituting these in (3.3) we find that

N 0 a * n z n satisfies (3.3) modulo O(z N +1
) for any N (in D ε and thus everywhere). Hence, the formal power series ∞ n=0 a * n z n satisfies (3.3), and thus a * n = a n for all n. In other words, (3.5) holds for z ∈ D ε . We have thus proved (vii) for a smaller domain D ε . Moreover, (vi) and (ix) follow.

Consider now instead a domain D ε := {z : |z| < ρ, | arg(z)| < π/2 -ε} in the right half-plane (with 0 < ε < π/2), and any solution A(z) there. We use (3.7), choosing z 0 with 0 < z 0 < ρ and integrating along the straight line from 0 to 1/z -1/z 0 . It is easily seen that, if |z| is small enough, then |w| ≤ z 0 < ρ and w = z/(1 -zt) ∈ D ε for all t on this line. We split the integral in (3.7) into two parts: one with |t| ≤ 1/(2|z|) and one with |t| > 1/(2|z|). In the second part, |e -t | ≤ e -sin(ε)/(2|z|) , while

|(1 -zt) γ | = |(z/w) γ | = O(|z| -γ + 1
) and the other factors are bounded; thus this part is O(z N +1 ) for any N . In the first part, |1 -zt| ≥ 1/2. We first substitute (for small |z| and with

F = R 3 e -G as above, F (z/(1 -zt)) = N m=0 c m (z/(1 -zt)) m + O(z N +1
), and then make a similar Taylor expansion of (1 -zt) γ-m . We combine this with a Taylor expansion of e G(z) and note that

1/(2z) 0 t k e -t dt = k! + O(z N +1
) for any k, and conclude that, for any choice of C in (3.7), A(z) has an asymptotic expansion A(z) ∼ ∞ n=0 a * n z n as z → 0 in D ε for some coefficients a * n . It follows by Cauchy's estimates, see Remark 2.1, that further A (z) ∼ ∞ n=0 a * n nz n-1 as z → 0 in D ε , for every ε < ε; consequently, arguing as for A 0 above, the formal power series ∞ n=0 a * n z n satisfies (3.3) and thus a * n = a n for all n. In other words, (3.5) holds for z ∈ D ε , for any solution A(z) of (3.3) in D ε . Differentiating (3.5) N times, we find, by Cauchy's estimates again, that D N A(z) → a N N ! as z → 0 in any strictly smaller domain D ε . Since ε ∈ (0, π/2) is arbitrary, this implies that (viii) holds. In particular, considering real z only, (v) follows.

We can now complete the proof of (vii). Fix ε < π/4 and let D + := {z ∈ D ε : z > 0} and D -:= {z ∈ D ε : z < 0}. The solution A 0 (z) is defined in D ε , and thus in D 1 := D ε ∩ D + . Choose any z 0 ∈ D 1 . Then A 0 (z) is given by (3.7) in D 1 , for some C; hence (3.7) defines a solution

A + (z) in D ε such that A + (z) = A 0 (z) in D ε ∩ D + . Similarly, there is a solution A -(z) in D ε such that A -(z) = A 0 (z) in D ε ∩ D -. (But note that A -= A + .) We extend A 0 to {z : |z| < δ, z / ∈ [0, ∞)} by defining A 0 (z) = A + (z) for z ∈ D + and A 0 (z) = A -(z) for z ∈ D -.
The claims in (vii) follow by the result for D ε proved above together with (viii) applied to D + and D -.

If we combine the solution A 0 on I -with any solution A(z) on I + , we can define functions f n on I by

f n (x) =      D n A 0 (x), x ∈ I -, a n n!, x = 0, D n A(x), x ∈ I + ;
by (v) and (vi), each function f n is continuous on I. Then f n (x) = f n (0) + x 0 f n+1 (y) dy, and thus each f n is differentiable with f n = f n+1 . Consequently, A := f 0 ∈ C ∞ (I). Further, A satisfies (3.3) on I \ {0}, and thus by continuity at 0 too. The claims (i)-(iv) now follow, using (v) and (vi).

We call A 0 the principal solution of (3.3).

Remark 3.3. It follows from the proof that we can define A 0 to be analytic in a disc {z : |z| < δ} cut along any given ray in the right half-plane. In particular, we can choose a ray in the fourth quadrant, say, and then A 0 is defined on (-δ, 0) ∪ (0, δ); if we further define A 0 (0) = a 0 , the restriction of this A 0 to (-δ, δ) is a C ∞ solution on an interval I as in (i). Moreover, we can regard A 0 (z) as a multivalued function; more precisely, we can regard A 0 (re iθ ) as a function of (r, θ) ∈ (0, ∞) × (-∞, ∞), and then (3.5) holds as r = |z| → 0 for -1 2 π + ε < θ < 5 2 π -ε, for any ε > 0.

The exponential generating function and Borel summation

We continue to study the recurrence (3.1), with α 1 > 0 as in Section 3, but consider instead of A(z) the exponential generating function

E(z) := ∞ n=0 a n z n n! . ( 4
* k = β k + Kα k , D K E(z) = K k=1 α k zD K-k+1 E(z) + β * k D K-k E(z) (4.2) or (1 -α 1 z)D K E(z) = K-1 j=1 (α K+1-j z + β * K-j )D j E(z) + β * K E(z). (4.3)
We use the standard method of introducing the vector-valued function

F (z) = (F k (z)) K-1 k=0 := (D k E(z)) K-1 k=0
and see that (4.3) is equivalent to

DF K-1 (z) = K-1 j=1 α K+1-j z + β * K-j 1 -α 1 z F j (z) + β * K 1 -α 1 z F 0 (z), (4.4) DF k (z) = F k+1 (z), 0 ≤ k ≤ K -2. (4.5) (When K = 1, we simply have F (z) = F 0 (z) = E(z), (4.5
) is vacuous and (4.4) is the same as (4.3).)

The system (4.4)-( 4.5) has a unique analytic solution with F (0) = (a k ) K-1 k=0 in any simply connected domain containing 0 but not α -1 1 , for example in

D E := C \ [α -1 1 , ∞); hence F 0 (z) is a solution to (4.3) in D E .
We let in the sequel E(z) denote this analytic extension of (4.1) to D E .

Remark 4.1. The example a n+1 = (n+β)a n , a 0 = 1, with E(z) = (1-z) -β , shows that in exceptional cases, E(z) might be analytic also at α -1 1 (and then E(z) is an entire function), and that it is also possible that E(z) has a pole at α -1 1 ; however, the same example shows that typically E(z) is singular at α -1 1 and that it is not possible to extend E(z) around α -1 1 as a single-valued function. See also Remark A.4.

The (generalized) Borel summation method for a (possibly divergent) series n a n consists of the following three steps, see [START_REF] Borel | Mémoire sur les séries divergentes[END_REF][START_REF] Borel | Lec ¸ons sur les séries divergentes[END_REF] and Hardy (1949, §8.5, 8.11). (Step (ii) is Hardy's generalization (Hardy, 1949, §8.11); Borel assumes that E(z) is entire.) (i) Define the exponential generating function E(z) by (4.1) (assuming it to exist at least for small |z|).

(ii) Extend (if necessary) E(z) analytically along the entire positive real axis (assuming this to be possible).

(iii) Define the Borel sum by

(B * ) ∞ n=0 a n := ∞ 0 E(x)e -x dx (4.6)
(assuming the integral to exist).

Theorem 4.2. Let (a n ) ∞ 0 be given by (3.1). There exists

δ > 0 such that if z ∈ D δ := {z ∈ C \ [0, ∞) : |z| < δ}, then the Borel sum (B * ) ∞ n=0 a n z n = ∞ 0 E(xz)e -x dx (4.7)
exists; furthermore, this sum is the principal solution of (3.3) defined in Section 3. Consequently, the Borel sum (B * ) ∞ n=0 a n z n solves (3.3) in D δ , and has there the asymptotic expansion ∼ ∞ n=0 a n z n and all derivatives are continuous also at 0, with the limit values a n there.

If ), and is an analytic solution of (3.3) there.

α k = 0 for 2 ≤ k ≤ K (in particular if K = 1), then we may take δ = ∞; i.e., then (B * ) ∞ n=0 a n z n exists in C \ [0, ∞
See Hardy (1949, Theorem 136) for a related general theorem on generalized Borel summability of asymptotic series. To prove Theorem 4.2 we begin with a simple estimate. Lemma 4.3. E(z) is an analytic function in D E , and satisfies the following estimates:

(i) |E(z)| ≤ e C|z| for some constant C and all z ∈ D E with |z| ≥ 2α -1 1 . (ii) If α k = 0 for 2 ≤ k ≤ K (in particular if K = 1), then |E(z)| ≤ e o(|z|) as z → ∞ with z ∈ D E .
Proof: We have already shown that E(z) extends to D E .

(i): The system (4.4)-(4.5) can be written

DF (z) = M (z)F (z) for a matrix M (z) with M (z) ≤ C for |z| ≥ 2α -1
1 and some constant C. Hence, for r ≥ r 0 := 2α -1 1 and 0 < θ < 2π, by Gronwall's lemma, see e.g. Revuz and Yor (1999, Appendix §1),

F (re iθ ) ≤ e C(r-r0) F (r 0 e iθ ) ,
and (i) follows (possibly with a larger C) since θ → F (r 0 e iθ ) can be extended to a continuous function on the closed interval [0, 2π], and thus is bounded.

(ii): Let ε > 0 and define

F ε (z) := (ε -k D k E(z)) K-1 k=0 .
We have, in analogy with (4.4)-(4.5) and using the assumption

α k = 0 for k ≥ 2, DF ε K-1 (z) = K-1 j=0 β * K-j 1 -α 1 z ε j+1-K F ε j (z), DF ε k (z) = εF ε k+1 (z), 0 ≤ k ≤ K -2, which can be written DF ε (z) = M ε (z)F ε (z), where M ε (z) ≤ 2ε for |z| large, say |z| ≥ R(ε). Consequently, by Gronwall's lemma again, |E(z)| ≤ |F ε (z)| ≤ Ce 2ε|z| for |z| ≥ R(ε) and some C = C(ε). Proof Proof of Theorem 4.2: Let z ∈ C \ [0, ∞).
The exponential generating function of the sequence (a n z n ) ∞ 0 is x → E(zx), which by Lemma 4.3 is analytic on [0, ∞) and there satisfies |E(zx)| ≤ C 1 e C2|z|x . Consequently, the integral (4.6) defining the Borel sum (B * ) ∞ n=0 a n z n becomes (4.7) and it exists at least for |z| < δ := C -1 2 ; we thus have

A(z) := (B * ) ∞ n=0 a n z n = ∞ 0 E(xz)e -x dx, z ∈ D δ . (4.8)
If α k = 0 for k ≥ 2, we may by Lemma 4.3(ii) choose C 2 arbitrarily small, and thus we may take

δ = ∞. Since E(z) is analytic in a disc {z : |z| < 2r} (with r := α -1 1 /2), we have, for any N , E(zx) = N n=0 a n x n z n /n! + O(|xz| N +1
) when |xz| ≤ r, and thus, for z ∈ D δ/2 ,

A(z) = r/|z| 0 N n=0 a n x n z n n! + O(|xz| N +1 ) e -x dx + ∞ r/|z| O e C2|z|x e -x dx = N n=0 ∞ 0 a n z n n! x n e -x dx + O(z N +1 ) + ∞ r/|z| O e x/2-x dx = N n=0 a n z n + O(z N +1 ), (4.9)
which shows the asymptotic expansion. A(z) is analytic in D δ , by (4.8) and standard arguments based on dominated convergence, and

A (z) = ∞ 0 xE (xz)e -x dx.
Consequently, using integration by parts K -k times in the first integral, and defining α 0 := 0, β0 :=

β * 0 := -1 for convenience, K k=1 α k z k+1 A (z) + K k=1 βk z k A(z) -A(z) = K k=0 ∞ 0 α k z k+1 xE (xz) + βk z k E(xz) e -x dx = K k=0 ∞ 0 D K-k x α k z k+1 xE (xz) + βk z k E(xz) e -x dx + K k=0 K-k-1 j=0 D j x α k z k+1 xE (xz) + βk z k E(xz) x=0
The integrated part is a polynomial r 1 (z) of degree at most K -1. Further,

D K-k x α k z k+1 xE (xz) + βk z k E(xz) = α k z K+1 xD K-k+1 E(xz) + β * k z K D K-k E(xz),
and summing we obtain

z K K k=0 α k xzD K-k+1 E(xz) + β * k D K-k E(xz) = 0
by (4.2); hence the last sum of integrals vanishes, and

K k=1 α k z k+1 A (z) + K k=1 βk z k A(z) -A(z) = r 1 (z) (4.10)
Consequently, A(z) satisfies the differential equation (4.2) in D δ , except that the polynomial r(z) is replaced by -r 1 (z). However, this means only that we have the differential equation for same recurrence (3.1) with some other initial values; however, the initial values are (by Theorem 3.2) given by the first terms of the asymptotic expansion, so by (4.9) they are just a 0 , . . . , a K-1 , and thus r(z) = -r 1 (z) and (4.10) equals (3.3). Hence, A(z) solves (3.3); since it is, by (4.9), a solution bounded as z 0 along the negative real axis, it is by Theorem 3.2 the principal solution A 0 (z).

5 Solving (1.1)

We return to the special case (1.1), for definiteness with a 0 = 0 and a 1 = 1 as in Section 1. The general solution A(z) to (1.3) is given in (1.6). Assume for simplicity α > 0. The general theory in Section 3 then shows that there is a principal solution A 0 (z) analytic in C \ [0, ∞) that is given, at least for z < 0, by taking z 0 = 0 and C = 0 in (1.6). We now make the change of variable t = α -1 (z -1 -w -1 ), or w = z/(1 -αzt), and obtain

A 0 (z) = z ∞ 0 (1 -αzt) -ν-1 e -t-tz 2 /(1-αzt) dt;
(5.1) this integral is an analytic function of z ∈ C \ [0, ∞), and the formula is thus valid for all such z, cf. (3.8).

(Recall that (3.8) was derived under the assumption α 1 = 1. In fact, if α = 1, we have in the notation of Section 3 γ = -ν, G(z) = z and R 3 (z) = z, and (5.1) coincides with (3.8). For a general α > 0, one can similarly derive (5.1) from (3.8) and the reduction in Remark 3.1.)

Consider for example -1 < z < 0, and fix N > 0. By two Taylor expansions,

A 0 (z) = z ∞ 0 (1 -αzt) -ν-1 N k=0 1 k! -tz 2 1 -αzt k e -t dt + O(z 2N +3 ) = N k=0 z 2k+1 k! ∞ 0 (1 -αzt) -ν-k-1 (-t) k e -t dt + O(z 2N +3 ) = N k=0 z 2k+1 k! ∞ 0 N =0 (ν + k + 1) ! (αzt) (-t) k e -t dt + O(z N +1 ) = N k=0 N =0 (-1) k α (ν + k + 1) (k + )! k! ! z 2k+ +1 + O(z N +1 ).
Since A 0 (z) has the asymptotic expansion (3.5) as z 0, we obtain by identifying coefficients the following solution to (1.1). (The argument has assumed α > 0, but since obviously each a n is a polynomial in α and β, the solution is valid for all α and β.)

Theorem 5.1. The solution of the recurrence (1.1) with a 0 = 0, a 1 = 1 is given by, for any complex α and β, and with ν := β/α (provided α = 0)

a n = (n-1)/2 k=0 (-1) k (n -k -1)! k! (n -2k -1)! α n-2k-1 (ν + k + 1) n-2k-1 = (n-1)/2 k=0 (-1) k n -k -1 k n-2k-1 j=1 ((k + j)α + β).
Knowing this solution, it should be easy to verify it directly by induction; we leave this as an exercise.

Remark 5.2. By continuity, the last formula is valid also in the case α = 0, where we thus obtain

a n = (n-1)/2 k=0 (-1) k n -k -1 k β n-2k-1 . (5.2)
We have also, by the standard argument for linear recurrences with constant coefficients (provided β = ±2),

a n = β + β 2 -4 n -β -β 2 -4 n β 2 -4 2 -n . (5.3)
Hence, the two formulas in (5.2) and (5.3) must agree. If we expand (5.3) into a polynomial, using the binomial theorem thrice, and compare coefficients of β n-2k-1 , we can see that this identity is equivalent to the hypergeometric sum j n 2j+1

j k = 2 n-2k-1 n-k-1 k . -1 0 1 2 3 -2 -1 x 1 0 -1 -1 0 -1 0 1 2x 0 1 0 -1 -x -2x 2 + 1 1 0 1 0 -1 -2x 2 -1 x 1 0 -1 3 -2x 2x 2 -1 2x 1 0 4 -6x 2 + 1 6x 3 -4x 6x 2 -1 3x 1 5 -24x 3 + 6x 24x 4 -18x 2 + 1 24x 3 -6x 12x 2 -1 4x
Tab. 1: Ψ(µ, ν; x) for some integers µ and ν.

For example, continuing (A.12)-(A.13),

Ψ(ν + 2, ν; x) = (ν + 1)x, Ψ(ν + 3, ν; x) = (ν + 1)(ν + 2)x 2 -1, Ψ(ν + 4, ν; x) = (ν + 1)(ν + 2)(ν + 3)x 3 -2(ν + 2)x, Ψ(ν + 5, ν; x) = (ν + 1) 4 x 4 -3(ν + 2) 2 x 2 + 1, Ψ(ν + 6, ν; x) = (ν + 1) 5 x 5 -4(ν + 2) 3 x 3 + 3(ν + 3)x.
Thus, typically the degree of Ψ(n + ν, ν; x) is exactly n -1 (when n ≥ 1), but when ν is a negative integer it may be smaller. We give some more examples in Tables 1234, where µ and ν are integers in Tables 1-3 and halfintegers in Table 4. (Note that for integer ν, J -ν(z) = (-1) ν J ν (z) and Y -ν(z) = (-1) ν Y ν (z); hence, Φ(µ, -ν; z) = (-1) ν Φ(µ, ν; z) and Ψ(µ, -ν; x) = (-1) ν Ψ(µ, ν; x), and similarly with the argument -µ if µ is an integer. This explains the symmetry in the tables, and the 0's when µ + ν = 0; these are special for integer arguments.) These tables are thus parts of doubly infinite arrays that satisfy the recurrences (A.10)-(A.11) in both directions; note that these arrays of polynomials and integers are constructed from the transcendental Bessel functions. (In the half-integer case the Bessel functions can be expressed in the trigonometric functions and powers of z.) We give also in Table 5 some examples from [START_REF] Sloane | The On-Line Encyclopedia of Integer Sequences[END_REF] of such integer sequences.

Remark A.4. As n → ∞, for any fixed ν and z, cf. Abramowitz and Stegun (1972, (9.3.1)),

J n+ν (z) ∼ (z/2) n+ν /Γ(n + ν + 1) and Y n+ν (z) ∼ -π -1 (z/2) -n-ν Γ(n + ν); hence |J n+ν (z)| → 0 and |Y n+ν (z)| → ∞ rapidly. Thus, by Theorem A.1, the solution a n of (1.1) satisfies |a n | ∼ c|α| n |Γ(n + ν)| ∼ cn ν-1 |α| n n!
for some c > 0, except in the special case a 1 J ν (2/α) = a 0 J ν+1 (2/α), when a n = cJ n+ν (2/α) → 0 rapidly. Hence, as claimed in Section 1, typically the generating function A(z) diverges for all z = 0, while the exponential generating function has radius of convergence |α| -1 , as seen more generally in Remark 4.1.

- Tab. 4: Ψ(µ, ν; 2) = 1 2 Φ(µ, ν; 1) for some half-integers µ and ν.

2 -1 0 1 2 3 4 5 6 7 -2 0 -1 1 1 0 -1 -3 -11 -52 -301 -1 1 0 -1 0 1 2 5 18 85 492 0 -1 1 0 -1 -1 -1 -2 -7 -33 -191 1 -1 0 1 0 -1 -2 -5 -18 -85 -492 2 0 -1 1 1 0 -1 -3 -11 -52 -301 3 1 -2 1 2 1 0 -1 -4 -19 -110 4 3 -5 2 5 3 1 0 -1 -5 -29 5
αn + β a n by Theorem A.1 A053983 2n + 1 Ψ(n + 1/2, 1/2; 2) -Ψ(n + 1/2, 3/2; 2) A053984 2n + 1 Ψ(n + 1/2, 1/2; 2) A058797 n + 1 Ψ(n + 1, 1; 1) -Ψ(n + 1, 2; 1) = Ψ(n + 1, 0; 1) A058798 n + 1 Ψ(n + 1, 1; 1) A058799 n + 3 Ψ(n + 3, 2; 1) A093985 2n Ψ(n, 0; 2) A093986 2n Ψ(n, 0; 2) -Ψ(n, 1; 2) A106174 2n + 2 Ψ(n + 1, 1; 2) A121323 2n + 3 Ψ(n + 3/2, 3/2; 2) A121351 3n + 4 Ψ(n + 4/3, 4/3; 3) A121353 3n + 1 Ψ(n + 1/3, 1/3; 3) A121354 3n + 2 Ψ(n + 2/3, 2/3; 3) Tab. 5: Some integer sequences satisfying (1.1) from [START_REF] Sloane | The On-Line Encyclopedia of Integer Sequences[END_REF].

Remark A.5. The Bessel functions of the third kind H

(1) ν (z) and H

(2) ν (z) are given by J ν (z) ± iY ν (z) (Abramowitz and Stegun, 1972, (9.1.3-4)), and thus

H (1) µ (z) H (1) ν (z) H (2) µ (z) H (2) ν (z) = 1 i 1 -i J µ (z) J ν (z) Y µ (z) Y ν (z)
.

Consequently, taking determinants and noting 1 i 1 -i = -2i, also Φ(µ, ν; z) = πi 2

H

(1) µ (z) H

(1) ν (z) H Remark A.6. The function J -ν satisfies the same Bessel equation z 2 f (z) + zf (z) + (z 2 -ν 2 )f (z) = 0 as J ν and Y ν , and it can be expressed in them by J -ν (z) = cos(νπ)J ν (z) -sin(νπ)Y ν (z) (Abramowitz and Stegun, 1972, (9.1.2)). It follows that, if µ -ν = n ∈ Z, J ν+n (z) J ν (z) J -ν-n (z) (-1) n J -ν (z) = (-1) n+1 sin(νπ) π Φ(ν + n, ν; z).

A more general formula expressing Φ(µ, ν; z) in Bessel functions J only is given in (B.3) below.

Remark A.7. The related recurrence a n+1 = (αn + β)a n + a n-1 can be solved similarly using the Bessel functions I ν and K ν ; alternatively, we may note that i n a n satisfies (1.1) (with iα and iβ instead of α and β), and thus the solution follows from Theorems 5.1 and A.1. In particular, it is easily seen that if a 0 = 0 and a 1 = 1, then the formulas in Theorem 5.1 hold if we omit the factor (-1) k . Examples of such sequences are A000806, A001053, A001517 and A001518 in [START_REF] Sloane | The On-Line Encyclopedia of Integer Sequences[END_REF], and the Bessel polynomials y n (x) that satisfy this recurrence with α = 2x, β = x (and thus ν = 1/2), and y 0 (x) = 1, y 1 (x) = x + 1, see [START_REF] Krall | A new class of orthogonal polynomials: The Bessel polynomials[END_REF] and [START_REF] Grosswald | Bessel Polynomials[END_REF].

Remark A.8. We do not know any closed form for the exponential generating function E(z) in general, but we observe that in the special case a n = J n (y) for some y = 0, when (1.1) holds with α = 2/y, β = 0 by (A.1),

E(z) := ∞ n=0 J n (y) n! z n = J 0 y 2 -2yz .
(This is an entire function, since J 0 is entire and even.) This can be verified by inserting the right-hand side into (4.3), using the Bessel equation. It would be interesting to know whether there are other cases when E(z) can be expressed using Bessel functions.

B Hypergeometric series

We have shown above (Theorem A.2) that if µ -ν ∈ Z, then Ψ(µ, ν; x) is a polynomial in x and, equivalently, Φ(µ, ν; z) is a polynomial in z -1 . In this appendix we will show this in another way by deriving the general formula (B.4) for Φ(µ, ν; z) in terms of hypergeometric functions. In the case µ-ν ∈ Z we show that this reduces to a polynomial in z -1 and we will recover Theorem A.2. However, the formula also shows that for general µ and ν no similar simplification is possible.

We begin with a product formula for Bessel functions (Abramowitz and Stegun, 1972, (9.1.14))

J µ (z)J ν (z) = z 2 µ+ν ∞ k=0 (µ + ν + k + 1) k Γ(µ + k + 1)Γ(ν + k + 1)k! - z 2 4 k . (B.1)
For µ, ν, µ + ν / ∈ Z <0 , this can also be written in terms of a hypergeometric function:

J µ (z)J ν (z) = (z/2) µ+ν Γ(µ + 1)Γ(ν + 1) 2 F 3 µ + ν + 1 2 , µ + ν + 2 2 ; µ + 1, ν + 1, µ + ν + 1; -z 2 . (B.2)
Remark B.1. Since, for ν / ∈ Z <0 , J ν (z) = (z/2) ν Γ(ν + 1) 0 F 1 (; ν + 1; -z 2 /4), this product formula is equivalent to the product formula for hypergeometric functions, for any α, β with α, β, α + β -1 / ∈ Z ≤0 : 0 F 1 (; α; x) 0 F 1 (; β; x) = 2 F 3 α + β -1 2 , α + β 2 ; α, β, α + β -1; 4x .

If ν / ∈ Z, then (Abramowitz and Stegun, 1972, (9.1.2))

Y ν (z) = cos(νπ)J ν (z) -J -ν (z) sin(νπ)

and thus the definition (A.3) yields, for µ, ν / ∈ Z, Φ(µ, ν; z) = -π sin(νπ) J µ (z)J -ν (z) + π sin(µπ) J ν (z)J -µ (z) + π cot(νπ) -cot(µπ) J µ (z)J ν (z). (B.3)

  Theorem 3.2. With conditions and notations as above, let I ⊆ R be a finite or infinite open interval such that 0 ∈ I and p(z) = 0 on I. Let I + := I ∩ (0, ∞) and I -:= I ∩ (-∞, 0). Then (i) There exists a differentiable function A(z) on I that solves (3.3) there. Any such solution is C ∞ on I and analytic on I \ {0} = I -∪ I + , but in general no solution is analytic at 0.(ii) The solution in (i) is unique on I -but not on I + .(iii) Any solution A(z) of (3.3) on I satisfies D n A(0)/n! = a n , n ≥ 0.

  .1) It is easily seen from (3.1), by induction, that |a n | = O(C n n!) for some C, and thus the sum in (4.1) converges and E(z) is an analytic function in some disc {z : |z| < R}. Differentiation of (4.1) yields D k E(z) = n a n+k z n /n!, and it follows that zD k+1 E(z) = n a n+k nz n /n!. Hence the recurrence relation (3.1) is equivalent to the differential equation, with β

1365-8050 c 2010 Discrete Mathematics and Theoretical Computer Science (DMTCS), Nancy, France

A Bessel functions

We use the standard Bessel functions of the first and second kinds J ν (z) and Y ν (z), see e.g. Abramowitz and Stegun (1972, Chapter 9) or Lebedev (1972, Chapter 5). Recall that these are analytic functions of z and ν for z ∈ C \ (-∞, 0] (the complex plane cut along the negative real axis) and ν ∈ C. (Alternatively, one can think of them as analytic functions of z in the covering space of C \ {0}, i.e., as analytic functions of log z ∈ C \ {0}, or as functions of (r, θ) ∈ (0, ∞) × (-∞, ∞).) The same applies to the functions Φ and Ψ defined below, so we assume throughout the appendices that z is in this domain (or in this covering space). (When µ is an integer, J µ (z) is an entire function of z, and there is no problem at all.)

The Bessel functions satisfy the recurrence relations

Hence, for a fixed z, the functions J ν (z) and Y ν (z), as functions of ν (restricted to a sequence ν 0 + n, n ∈ N), satisfy a recurrence relation of the type studied in this paper. Conversely, the solution to (1.1) can easily be expressed in these Bessel functions by fitting the constants and initial values. We find it convenient to define and Stegun, 1972, (9.1.27)), the formula for the Wronskian W {J ν , Y ν } yields, as stated in Abramowitz and Stegun (1972, (9.1.16)),

It follows immediately from (A.1) and (A.2) that

Hence, for fixed z, the function Φ is a function of two (complex) variables that satisfy the type of recurrence studied here in each variable. To write the recurrence in a more convenient form, we define, for (complex

Further, Ψ(µ, ν; x) = -Ψ(ν, µ; x) and

Consequently, for any fixed ν and x = 0, the sequence Ψ(ν + n, ν; x), n ∈ N, satisfies (1.1) with α = x and β = νx, and the initial conditions a 0 = 0, a 1 = 1 used in Section 1. We just as easily find the general solution:

Theorem A.1. The general solution to the recurrence (1.1),

In particular, the solution with a 0 = 0, a 1 = 1 is

Proof: These expressions satisfy the recurrence by (A.10) and the initial conditions by (A.12), (A.13) and

Combining the two solutions of the same recurrence found in Theorems 5.1 and A.1, we find the following formulas.

Theorem A.2. For any complex ν and n ∈ N,

a polynomial of degree at most n -1 in x, and

, a polynomial of degree at most n in z -1 .

Remark A.3. These polynomials are called Lommel polynomials, see [START_REF] Lommel | Zur Theorie der Bessel'schen Functionen[END_REF] and Watson (1922, §9•6). In the notation of [START_REF] Watson | A Treatise on the Theory of Bessel Functions[END_REF], Ψ(n + ν, ν; x) = R n-1,ν+1 (2/x). For combinatorial properties of Lommel polynomials, see also Flajolet and[START_REF] Flajolet | Nonoverlapping partitions, continued fractions, Bessel functions and a divergent series[END_REF][START_REF] Feinsilver | Combinatorial interpretation and operator calculus of Lommel polynomials[END_REF][START_REF] Feinsilver | Combinatorial interpretation and operator calculus of Lommel polynomials[END_REF].

Svante Janson

The product formula (B.2) and Γ(y)Γ(1 -y) = π/ sin(yπ) yield, if µ, ν, µ ± ν / ∈ Z, the formula for Φ (and thus for Ψ) promised above:

In general, no simplification seems possible, since the three terms on the right hand side have expansions with powers of z that are congruent (mod 2Z) to, respectively, µ -ν, ν -µ and µ + ν. In the special (limiting) case µ -ν ∈ Z, however, cot(νπ) = cot(µπ), so the third term disappears. Moreover, if, say, µ -ν = n ∈ Z ≥0 , but still µ, ν / ∈ Z, then the result can be written (simplest seen from (B.3) and (B.1), noting that sin(µπ) = (-1) n sin(νπ)),

It is easy to see that the term with k = m, say, in the first sum is cancelled by the term with k = m + n in the second. Hence

, the polynomial of degree n in z -1 found in Theorem A.2 by a different method using our generating functions. We have thus found another proof of Theorem A.2. (The argument above assumes ν / ∈ Z, but since the final result is a polynomial in ν, and the left hand side is analytic in ν, the result holds for all ν by continuity.)