# Automatic classification of restricted lattice walks.

1 Introduction

There is a strange phenomenon about the generating functions that count lattice walks restricted to the quarter plane: depending on the choice of the set [??] [subset equal to] {[??], [left arrow], [??], [up arrow], [??], [right arrow], [??], [down arrow]} of admissible steps, the generating function is sometimes rational, sometimes algebraic [but not rational], sometimes D-finite [but not algebraic], and sometimes not even D-finite. This is quite in contrast to the corresponding problem in 1D, where the generating functions invariably are algebraic [3]. Much progress was made recently on understanding why this is so, and only very recently, Bousquet-Melou and Mishna [9] have announced a classification of all the 256 possible step sets into algebraic, transcendental D-finite, and non-D-finite cases, together with proofs for the algebraic and D-finite cases and strong evidence supporting the conjectured non-D-finiteness of the others.

As usual, a power series S(t) [member of] Q[[t]] is called algebraic if there exists a bivariate polynomial P(T, t) in Q[T, t] such that P(S(t), t) = 0, and transcendental otherwise. Also as usual, a power series S(t) is called D-finite if it satisfies a linear differential equation with polynomial coefficients. (Every algebraic power series is D-finite, but not vice versa.) At first glance, it might seem easy to prove that a power series is algebraic or D-finite: just come up with an appropriate equation, and then verify that the series satisfies this equation. But as far as lattice walks are concerned, most proofs given so far are indirect in that they avoid exhibiting the equation explicitly but merely are satisfied showing its existence. This is probably so because the equations appearing in this context are often too big to be dealt with by hand.

Nevertheless, it is interesting to know the equations explicitly, because they provide a standard canonical representation for a series, from which lots of further information can be extracted in a straightforward manner. By applying a well-known technique from computer algebra (in modern fashion, cf. Section 2), we have systematically searched for differential equations and algebraic equations that the series counting the walks in the quarter plane satisfy. These are given in Section 3. We have also made a first step towards classifying walks in [Z.sup.3] confined to the first octant (cf. Section 4) by considering all step sets [??] with up to five elements, and performed a systematic search for equations of the corresponding series. More than 2000 hours of computation time have been spent in order to analyze about 3500 different sequences.

We do not provide proofs that the equations we found are indeed correct, but the computational evidence in favor of our equations is striking. We have no doubt that all the equations we found are correct. In principle, it would be possible to supplement the "automatically guessed" equations by computer proofs in a systematic fashion, using techniques that have recently been applied to some special cases [25, 24, 6]. But we found that the computational cost for performing these automated proofs would be by far higher than what was needed for the mere discovery.

2 Methodology

To study generating functions for lattice walks, we follow a classical scheme in experimental mathematics. It is based on the following steps: (S1) computation of high order expansions of generating power series; (S2) guessing differential and/or algebraic equations satisfied by those power series; (S3) empirical certification of the guessed equations (sieving by inspection of their analytic, algebraic and arithmetic properties); (S4) rigorous proof, based on (exact) polynomial computations.

In what follows, we only explain Steps (S1), (S2) and (S3). A full description of Step (S4) is given in [6]. By way of illustration, we choose an example requiring computations with human-sized outputs, namely the classical case, initially considered by Kreweras [27, 7, 8], of walks in the quarter plane restricted to the step set [??] = {[left arrow], [??], [down arrow]}.

2.1 Basic Definitions and Facts

We focus on 2D and 3D lattice walks. The 2D walks that we consider are confined to the quarter plane [N.sup.2], they join the origin of [N.sup.2] to an arbitrary point (i, j) e [N.sup.2], and are restricted to a fixed subset [??] of the step set {[??], [left arrow], [??], [up arrow], [??], [right arrow], [??], [down arrow]}. If f (n; i, j) denotes the number of such walks of length n (i.e., using n steps chosen from [??]), the sequence f (n; i, j) satisfies the multivariate recurrence with constant coefficients

f (n + 1; i, j)= [summation over ((h, k)[member of][??])] f(n; i - h, j - k) for n, i, j [greater than or equal to] 0. (1)

Together with the appropriate boundary conditions

f(0;0,0) = 1 and f(n; i, j) = 0 if i < 0 or j < 0 or n < 0,

the recurrence relation (1) uniquely determines the sequence f(n; i, j). As is customary in combinatorics, we let

F(t; x, y) = [[summation].sub.n[greater than or equal to]0]( [[summation].sub.i, j[greater than or equal to]0]f(n; i, j)[x.sup.i][y.sup.i])[t.sup.n]

be the trivariate generating power series of the sequence f (n; i, j). As f (n; i, j) = 0 as soon as i > n or j > n, the inner sum is actually finite, and so we may regard F(t; x, y) as a formal power series in t with polynomial coefficients in Q[x, y].

Specializing F(t; x, y) to selected values of x and y leads to various combinatorial interpretations. Setting x = y = 1 yields the power series F(t; 1,1) whose coefficients count the total number of walks with prescribed number of steps (and arbitrary endpoint); the choice x = y = 0 gives the series F(t; 0, 0) whose coefficients count the number of walks returning to the origin; setting x = 1, y = 0 yields the power series whose coefficients count the number of walks ending somewhere on the horizontal axis, etc.

By [10, Th. 7], multivariate sequences that satisfy recurrences with constant coefficients have moderate growth, and thus their generating series are analytic at the origin. The next theorem refines this result in our context.

Theorem 1 The following inequality holds

f(n; i, j) [less than or equal to] [[absolute value of [??]].sup.n] for all (i, j, n) [member of] [N.sup.3]. (2)

In particular, the power series F(t; 0, 0), F(t; 1,0), F(t; 0, 1) and F(t; 1, 1) are convergent in C[[t]] at t = 0 and their radius of convergence is at least 1/[absolute value of [??]].

Proof: The total number of unrestricted n-step walks starting from the origin is [absolute value of [[??].sub.n]], so the number of walks restricted to a certain region is bounded by this quantity. This implies that the coefficient of [t.sup.n] in F(t; 1,1) is at most [[??].sub.n]]. The bound also applies to the coefficient of [t.sup.n] in F(t; [alpha], [beta]) for [alpha], [beta] [member of] {0, 1}, as these series count walks which are subject to further restrictions.

2.1.1 D-finite generating series of walks are G-functions

A power series S(t) = [[summation].sub.n[greater than or equal to]0] [a.sub.n][t.sup.n] in Q[[t]] is called a G- function (i) if (a) it is D-finite; (b) its radius of convergence in C[[t]] is positive; (c) there exists a constant C > 0 such that for all n [member of] N, the common denominator of [a.sub.0], ..., [a.sub.n] is bounded by [C.sup.n].

Examples of G-functions are the power series expansions at the origin of log(1 -1) and [(1 - t).sup.[alpha]] for [alpha] [member of] Q. More generally, the Gauss hypergeometric series [sub.2][F.sub.1]([alpha], [beta], [gamma]; t) with rational parameters [alpha], [beta], [gamma], is also a G-series [17]. A celebrated theorem of Eisenstein assures that any algebraic power series must be a G-function (if S is algebraic, there exists an integer C [member of] N such that [a.sub.n][C.sup.n+1] is an integer for all n.) The fact that G-functions arise frequently in combinatorics was recently pointed out by Garoufalidis [20].

G-functions enjoy many remarkable properties. Chudnovsky [14] proved that the minimal order differential equation satisfied by a G-series must be globally nilpotent (see Section 2.4.4 below for the definition and an algorithmic use of this notion). By a theorem of Katz and Honda [22, 21], the global nilpotence of a differential operator implies that all of its singular points are regular singular points with rational exponents. See also [1, 13, 17] for more details on this topic.

Theorem 2 Let S(t) be one of the power series F (t; 0, 0), F (t; 1, 0), F (t; 0, 1) and F (t; 1, 1). If S is D-finite, then S is a G-series. In particular, its minimal order homogeneous linear differential equation is Fuchsian and it has only rational exponents. Moreover, the coefficient sequence of S(t) is asymptotically equivalent to a sum of terms of the form [kappa][[rho].sup.n][n.sup.[alpha]][(log n).sup.[beta]] for some constants [kappa] [member of] R, [alpha] [member of] Q, [rho] [member of] [bar.Q], and [beta] [member of] N.

Proof: The conditions (a) and (c) in the definition of a G-function are clearly satisfied. The only non-trivial point is the fact that the series S has a positive radius of convergence in C. This follows from Theorem 1. The Fuchsianity of the minimal equation for S, and the rationality of its exponents, follow by combining the results by Katz, Honda and Chudnovsky cited above. The claim on the asymptotics of the coefficients of S(t) is a consequence of [20, Prop. 2.5].

For 3D walks, the definitions are analogous. The trivariate power series F(t; x, y) is simply replaced by the generating series G(t; x, y, z) [member of] Q[x, y, z][[t]] of the sequence g(n; i, j, k) that counts walks in [N.sup.3] starting at (0, 0, 0) and ending at (i, j, k) [member of] [N.sup.3]. Note that the appropriate versions of Theorems 1 and 2 hold; in particular, the generating series of octant walks G(t; 1, 1, 1) is a G-series whenever it is D-finite.

2.2 Computing large series expansions

The recurrence (1) can be used to determine the value of f (n; i, j) for specific integers n, i, j [member of] N. Theorem 1 implies that f (n; i, j) is a non-negative integer whose bit size is at most O(n). If N [member of] N, the values f (n; i, j) for 0 [less than or equal to] n, i, j [less than or equal to] N can thus be computed altogether by a straightforward algorithm that uses O([N.sup.3]) arithmetic operations and [??]([N.sup.4]) bit operations. (We assume that two integers of bit-size N can be multiplied in [??](N) bit operations; here, the soft-O notation [??]() hides logarithmic factors.) The memory storage requirement is proportional to [N.sup.3]. The same is also true for the truncated power series [F.sub.N] = F(t; x, y) mod [t.sup.N]. For our experiments in 2D, we have chosen N = 1000. With this choice, the computation of the f (n; i, j) is the step which consumes by far the most computation time in our calculations. (ii)

Example 1 The Kreweras walks satisfy the recurrence

f (n + 1, i, j) = f (n, i + 1, j) + f (n, i, j + 1) + f (n, i - 1, j - 1) for n, i, j [greater than or equal to] 0,

which allows the computation of the first terms of the series F (t; x, y)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

and also the first terms of the generating series F (t; 1,1) for the total number of Kreweras walks

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

In the 3D case, the values g(n; i, j, k) for 0 [less than or equal to] n, i, j, k [less than or equal to] N can be computed in O([N.sup.4]) arithmetic operations, [??]([N.sup.5]) bit operations and O([N.sup.4]) memory space. In practice, we found that computing G mod [t.sup.N] with N = 400 is feasible.

2.3 Guessing

Once the first terms of a power series are determined, our approach is to search systematically for candidates of linear differential equations or of algebraic equations which the series may possibly satisfy. This technique is classical in computer algebra and mathematical physics, see for example [11, 31, 28]. Differential and algebraic guessing procedures are available in some computer algebra systems like Maple and Mathematica.

2.3.1 Differential guessing

If the first N terms of a power series S [member of] Q[[t]] are available, one can search for a differential equation satisfied by S at precision N, that is, for an element L in the Weyl algebra Q[t]<[D.sub.t]> of differential operators in the derivation [D.sub.t] = d/dt with polynomial coefficients in t, such that

L(S) = [c.sub.r](t)[S.sup.(r)](t) + ... + [c.sub.1](t)S'(t) + [c.sub.0](t)S(t) = 0 mod [t.sup.N]. (3)

Here, the coefficients [c.sub.0](t), ..., [c.sub.r](t) [member of] Q[t] are not simultaneously zero, and their degrees are bounded by a prescribed integer d [greater than or equal to] 0. By a simple linear algebra argument, if d and r are chosen such that (d + 1)(r + 1) > N, then such a differential equation always exists. On the other side, if d, r and N are such that (d + 1)(r + 1) [much less than] N, the equation (3) translates into a highly over-determined linear system, so it has no reason to possess a non-trivial solution.

The idea is that if the given power series S(t) happens to be D-finite, then for a sufficiently large N, a differential equation of type (3) (thus satisfied a priori only at precision N) will provide a differential equation which is really satisfied by S(t) in Q[[t]] (i.e., at precision infinity). In other words, the D-finiteness of a power series can be (conjecturally) recognized using a finite amount of information.

Given the values d, r, N, and the first N terms of the series S, a candidate differential equation of type (3) for S can be computed by Gaussian elimination in O([N.sub.3]) arithmetic operations and [??]([N.sup.4]) bit operations. Actually, a modular approach is preferred to a direct Gaussian elimination over Q. Precisely, the linear algebra step is performed modulo several primes p, and the results (differential operators modulo p) are recombined over Q via rational reconstruction based on an effective version of the Chinese remainder theorem. (See [23] for an implementation of this technique in Mathematica.)

If no differential equation is found, this definitely rules out the possibility that a differential equation of order r and degree d exists. This does not, however, imply that the series at hand is not D-finite. It may still be that the series satisfies a differential equation of order higher than r or an equation with polynomial coefficients of degree exceeding d.

Asymptotically more efficient guessing algorithms exist, based on fast Hermite-Pade approximation [4] of the vector of (truncated) power series [S, S', ..., [S.sup.(r)]]; they have arithmetic complexity quadratic or even softly-linear in N. Such sophisticated algorithms were not needed to obtain the results of this paper, but they have provided crucial help in the treatment of examples of critical sizes (e.g. guessing with higher values of d, r, N and/or over a parametric base field like Q(x) instead of Q) needed for the proof in [6].

Example 2 (continued) N = 100 terms of the generating series F(t; 1, 1) of the total number of Kreweras walks are sufficient to conjecture that F(t; 1,1) is D-finite, since it verifies the differential equation [L.sub.,1,1](F(t; 1, 1)) = 0 mod [t.sup.N], where

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. (4)

Thus, with high probability, F(t; 1,1) verifies the differential equation [L.sub.1,1](F(t; 1, 1)) = 0.

Sometimes (see Section 2.4.4) one needs to guess the minimal-order differential equation [L.sub.min](S) = 0 satisfied by the given generating power series. Most of the time, the choice (d, r) of the target degree and order does not lead to this minimal operator. Worse, it may even happen that the number of initial terms N is not large enough to allow the recovery of [L.sub.min], while these N terms suffice to guess non-minimal order operators. (The explanation of why such a situation occurs systematically was given in [5], for the case of differential equations satisfied by algebraic functions.) A good heuristic is to compute several non-minimal operators and to take their greatest common right divisor; generically, the result is exactly [L.sub.min].

As a final general remark, let us point out that a power series satisfies a linear differential equation if and only if its coefficients satisfy a linear recurrence equation with polynomial coefficients. A recurrence equation can be computed either from a differential equation, or it can be guessed from scratch by proceeding analogously as described above for differential equations.

Example 3 (continued) N = 100 terms of the series S(t) = F(t; 1, 1) suffice to guess that its coefficients satisfy the order-6 recurrence

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

2.3.2 Algebraic guessing

If the first N terms of a power series S [member of] Q[[t]] are available, one can also search for an algebraic equation satisfied by S at precision N, that is, for a bivariate polynomial P(T, t) in Q[T, t] such that

P(S(t),t) = [c.sub.r](t)S[(t).sup.r] + ... + [c.sub.1](t)S(t) + [c.sub.0](t) = 0 mod [t.sup.N]. (5)

A similar discussion shows that candidate algebraic equations of type (5) for S can be "guessed" by performing either Gaussian elimination or Hermite-Pade approximation on the vector [1, S, ..., [S.sup.r]], followed by a gcd computation in Q[T, t] applied to two (or more) different guesses.

Example 4 (continued) N = 100 terms of the series S(t) = F (t; 1, 1) counting the total number of Kreweras walks suffice to guess that F (t; 1, 1) is very probably algebraic, namely solution of the bivariate polynomial

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. (6)

2.4 Empirical certification of guesses

Once discovered a differential equation (3) or an algebraic equation (5) that the power series S(t) seems to satisfy, we inspect several properties of these equations, in order to provide more convincing evidence that they are correct. These properties have various natures: some are computational features (moderate bit sizes), others are algebraic, analytic and even arithmetic properties. We check them systematically on all the candidates; if they are verified, as in the Kreweras example, this offers striking evidence that the guessed equations are not artefacts.

2.4.1 Size sieve: Reasonable bit size

The differential equation (3) has typically much lower bit size than a differential equation produced by the same guessing procedure applied to the same order, degree and precision, but to an arbitrary series having coefficients of bit-size comparable to that of S(t). A similar observation holds for the algebraic equation (5).

Example 5 (continued) If we perturb the coefficients of S(t) = F(t; 1, 1) by just adding a random integer between -100 and 100 to each of its coefficients, then the differential guessing procedures at order r = 4, degree d = 9 and precision N = 100 will either give no result (the over-determined system approach) or produce fake candidates (the Hermite-Pade approach) with polynomial coefficients in t, whose coefficients in Q have numerators and denominators of about 500 decimal digits each, instead of 4 digits for [L.sub.1,1].

2.4.2 Algebraic sieve: High order series matching

The equations (3) and (5) were obtained starting from N coefficients of the power series S(t). They are therefore satisfied a priori only modulo [t.sup.N]. We compute more terms of S(t), say 2N, and check whether the same equations still hold modulo [t.sup.2N]. If this is the case, chances increase that the guessed equations also hold at infinite precision.

2.4.3 Analytic sieve: Singularity analysis

By Theorem 2, the minimal order operators for power series like S(t) = F(t; 0, 0) and S(t) = F(t; 1, 1) must have only regular singularities (including the point at infinity) and their exponents must be rational numbers.

Example 6 (continued) The differential operator [L.sub.1,1] is Fuchsian. Indeed, a (fully automated) local singularity analysis shows that the set of its singular points {-1, 0, [infinity], 1/3, 4/3, -1/6 (1 [+ or -] i[square root of 3])} is formed solely of regular singularities. Moreover, the indicial polynomials of [L.sub.1,1] are, respectively: t(t - 1)(t - 2)(2t - 1), t(t - 1)(2t + 1)(t + 1), (t - 5)(t - 1)(t - 2)(t - 4), (t + 1)t(4t - 1)(4t + 1), t(t - 1)(t - 2)(t - 4), and t(t - 2)(2t - 3)(t - 1). Their roots are the rational exponents of the singularities.

2.4.4 Arithmetic sieve: G-series and global nilpotence

Last, but not least, we check an arithmetic property of the guessed differential equations by exploiting the fact that those expected to arise in our combinatorial context are very special.

Indeed, by a theorem due to the Chudnovsky brothers [14], the minimal order differential operator L [member of] Q[t]<[D.sub.t]> killing a G-series enjoys a remarkable arithmetic property: L is globally nilpotent. By definition, this means that for almost every prime number p (i.e., for all with finitely many exceptions), there exists an integer [mu] [greater than or equal to] 1 such that the remainder of the Euclidean (right) division of [D.sup.p[mu].sub.t] by L is congruent to zero modulo p [21, 16].

From a computational view-point, a fine feature is that the nilpotence modulo p is checkable. If r denotes the order of L, let [M.sub.p] be the p-curvature matrix of L, defined as the r x r matrix with entries in Q(t) whose (i, j) entry is the coefficient of [D.sup.j-1.sub.t] in the remainder of the Euclidean (right) division of [D.sup.p+i- 1.sub.t] by L. Then, L is nilpotent modulo p if and only if the matrix [M.sub.p] is nilpotent modulo p [16, 32].

In combination with Theorem 2, this yields a fast algorithmic filter: as soon as we guess a candidate differential equation satisfied by a generating series which is suspected to be a G-series (e.g. by F(t; 1, 1)), we check whether its p-curvature is nilpotent, say modulo the first 50 primes for which the reduced operator L mod p is well-defined. If the p-curvature matrix of L is nilpotent modulo p for all those primes p, then the guessed equation is, with very high probability, the correct one.

We push even further this arithmetic sieving. A famous conjecture, attributed to Grothendieck, asserts that the differential equation L(S) = 0 possesses a basis of algebraic solutions (over Q(x)) if and only if its p-curvature matrix [M.sub.p] is zero modulo p for almost all primes p. Even if the conjecture is, for the moment, fully proved only for order one operators and partially in the other cases [13], we freely use it as an oracle to detect whether a guessed differential equation has a basis of algebraic solutions. For instance, the computation of the p-curvature of an order 11 differential operator with polynomial coefficients of degree 96 in t, was one of the key points in our discovery [6] that the trivariate generating function for Gessel walks is algebraic.

Example 7 (continued) The 5-curvature matrix [M.sub.5](t) of the differential operator [L.sub.1,1] in (4) has the form 1/d(t)[M.sub.5](t), where d(t) = [(3t - 1).sup.7][t.sup.6][(t + 1).sup.5][(9[t.sup.2] + 3t + 1).sup.5](3t - 4) and [[??].sub.5](t) is a 4 x 4 matrix with polynomial entries in Q[t] of degree at most 27. The characteristic polynomial X[M.sub.5] of [M.sub.5] reads

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

where [N.sub.0], [N.sub.1], [N.sup.2], [N.sub.3] are irreducible polynomials in Z[t], of degree, respectively, 21, 26, 26, 21 and with coefficients having at most 20 decimal digits.

The polynomial X[M.sub.5] obviously equals [T.sup.4] modulo p = 5, so the 5-curvature of [L.sub.1,1] is nilpotent (but not zero (iii)) modulo 5. In fact, for all the primes 7 [less than or equal to] p < 100, the p-curvature matrix of [L.sub.1,1] is also nilpotent modulo p; it is even zero modulo p. Under the assumption that Grothendieck's conjecture is true, this indicates that [L.sub.1,1] admits a basis of algebraic solutions, and so provides independent evidence that also S(t) = F(t; 1, 1) is algebraic.

3 Empirical Results in 2D

In this section, we consider the total number of walks only, i.e., the generating function F(t; 1, 1). Because of symmetries, the 256 possible step sets give rise to 92 different sequences only. By inspection of the first N = 1000 terms, we found that 36 of them appear to be D-finite: 19 are algebraic and 17 are transcendental. The D-finite step sets, together with the sizes of the equations we discovered, are listed in Table 1 in the appendix. (There, and below, step sets are represented by compact pictograms, e.g. [??] for [??] = {[left arrow], [??], [down arrow]}.)

3.1 Combinatorial Observations

Our classification matches the results of Bousquet-Melou and Mishna [9]: for every sequence they prove D-finite our software found a recurrence and a differential equation, and whenever a series is algebraic indeed, our programs recognized it. Moreover, we found no recurrence or differential equation for any step set conjectured non-D-finite by Bousquet-Melou and Mishna. This strengthens the evidence in favor of the conjectured non-D-finiteness of these cases.

3.2 Algebraic Observations

All but two of the minimal polynomials of the algebraic series share the property that they define a curve of genus 0. As a consequence, there exists a rational parametrization in all these cases. For example, for the Kreweras step set [??], the minimal polynomial [P.sub.1,1] given in (6) defines a curve parameterized by

T(u) = ([u.sup.2] + 24u + 151)a(u)/(u + 9)([u.sup.2] + 24u + 147) and t(u) = 2/a(u),

where a(u) = ([u.sup.6] + 66[u.sup.5] + 1827[u.sup.4] + 27180[u.sup.3] + 229431[u.sup.2] + 1042866u + 1995717)/(u+11)[([u.sup.2]+22u+125).sup.2], i.e., for these rational functions we have

[P.sub.1,1](T(u), t(u)) = 0.

The two algebraic series that do not admit a rational parametrization belong to the step sets [??] (reverse Kreweras) and [??] (Gessel's). Their genus is 1.

Another feature of the series which we found to be algebraic is that they all admit closed forms in terms of (nested) radical expressions. For example, for the Kreweras step set, we find that F(t; 1, 1) is equal to

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

where i = [square root of -1] and [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Such representations can be found by appealing to the built-in equation solvers of Maple and Mathematica applied to the equation [P.sub.1,1] = 0. Both features are remarkable because, among all algebraic power series, only a few are rationally parameterizable or expressible in terms of radicals.

Also the transcendental D-finite series appear to have some special properties. Being D-finite, these series are annihilated by some linear differential operator

L = [c.sub.0](t) + [c.sub.1](t)[D.sub.t] + ... + [c.sub.r](t)[D.sup.r.sub.t] [member of] Q[t]<[D.sub.t]>).

According to the DFactor command from Maple's DEtools package, all the operators can be factorized into a product of one irreducible operator of order 2 and several operators of order 1. As all the operators are globally nilpotent, so are all their factors [16, 17].

We can therefore expect that every solution of these factors can be written as a sum of terms of the form

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], (7)

where R and Z are rational functions in Q(t) and [alpha], [beta], [gamma], [delta] are rational numbers. Indeed, Dwork [16, Item 7.4] has conjectured that any globally nilpotent second order differential equation has either algebraic solutions or is gauge equivalent to a weak pullback of a Gauss hypergeometric differential equation with rational parameters. This conjecture was disproved by Krammer [26] and recently by Dettweiler and Reiter [15]; the counter-examples given in these papers require involved tools in algebraic geometry (arithmetic triangle groups, systems associated to periods of Shimura curves, ...)

We are therefore in a win-win situation: either the second order operators appearing as factors of our operators admit only solutions which are indeed sums of terms of the form (7), or there is a simple combinatorial counter-example to Dwork's conjecture. Let us illustrate this on one of the most simple examples, the step set [??]. We find here the differential operator

4(32[t.sup.2] - 12t - 1) + 4(8t - 1)(20[t.sup.2] - 3t - 1)[D.sub.t] + t(4t - 1)(112[t.sup.2] - 5)[D.sup.2.sub.t] + [t.sup.2][(4t - 1).sup.2](4t + 1)[D.sup.3.sub.t]

which Maple factors into

(2(192[t.sup.3] - 56[t.sup.2] - 6t + 1) +4(24[t.sup.2] - 1)(4t - 1)t[D.sub.t] + [(4t - 1).sup.2](4t + 1)[t.sup.2][D.sup.2.sub.t])(1/t + [D.sub.t]).

With the help of Maple's built-in differential equation solver (the dsolve command), it can be found that the differential operator gives rise to the representation

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

(Incidentally, this solution can also be expressed in terms of elliptic functions.) We believe that all the transcendental D-finite generating functions for any step set admit a representation as (a nested integral of) such an expression. The solvers of Maple and Mathematica, however, are able to discover such a representation only in the simplest cases. (Note that at present, no complete algorithm is known that is capable of finding general pullback representations.)

3.3 Analytic Observations

By Theorem 2, all the coefficient sequences grow like [kappa][n.sup.[alpha]] [[rho].sup.n] log[(n).sup.[beta]] for some constants [kappa], [rho], [alpha], [beta] (we only care about the dominant part of their asymptotic expansions). From the differential equation or the recurrence equation, we can determine [rho], [alpha], and [beta] exactly as roots of characteristic polynomials and indicial equations, respectively. (See [33, 19] on how this is done.) We find that [beta] = 0 in all cases. Knowing the recurrence, we can also compute easily tens of thousands of sequence terms. With the help of convergence acceleration techniques [12] applied to so many terms, it is possible to determine the remaining constant k to an accuracy of thirty digits or more. With that many digits, it makes sense to search systematically for potential exact expressions of these constants using Plouffe's inverter [30] and/or algorithms like LLL and PSLQ [2]. We actually found "closed form" expressions for all these constants. They are included in Table 1 in the appendix.

By Theorem 1, the numbers [rho] are bounded by the cardinality of the step set [??]. It turns out that [rho] = [absolute value of [??]] unless the vector sum of the elements of the step set points outside the first quadrant. In these cases, [rho] is an algebraic number of degree 2 (e.g., p = 1 + 2[square root of 2] for the step set [??]). For [alpha], we found only non- positive numbers. Note that [alpha] being a negative integer implies that the corresponding series is transcendental [18].

All the constants [kappa] have the form [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], where the [[phi].sub.i] are usually small integers, the [e.sub.i] are rational numbers, and u is 1/[pi] if F(t; 1, 1) is transcendental, and 1/[GAMMA]([alpha] + 1) if F(t; 1, 1) is algebraic.

There are some cases where the [[phi].sub.i] are not integers. Among them, very strange is only the case of the step set [??], for which we found r = 1, [e.sub.o] = 7/2, [e.sub.1] = 1/2, [rho] = 2 + 2[square root of 6] and [[phi].sub.1] = (1137 + 468[square root of 6])/152000. This last number may look like a guessing artefact at first glance, but we trust in its correctness, because the number of correct digits exceeds by far the number of correct digits to be expected from an artefact.

4 Empirical Results in 3D

We have investigated walks in three dimensions confined to the first octant with step sets of up to five elements. A priori, there are 83682 such step sets, and they give rise to 3334 different sequences. Of those, we have computed the first N = 400 terms of the generating function G(t; 1,1,1) of general walks, and searched for potential differential equations, algebraic equations, and recurrence equations. We found that 134 sequences appear to be D-finite, and among those, 50 appear to be algebraic.

4.1 Combinatorial Observations

For some of the sequences, it can be realized that their D-finiteness or algebraicity is a consequence of the D-finiteness or algebraicity of a certain 2D walk. For example, the sequence corresponding to the step set

[??] 1, 4, 17, 75, 339, 1558, 7247, 34016, 160795, 764388, ... (A026378)

is readily seen to be D-finite, since it may be regarded as a variation of the 2D step set [??] in which the step [up arrow] appears in two copies and empty steps are allowed. (Here and below, a three dimensional step set is depicted in three separate slices: first the arrows tops of the forms (x, y, -1), then (x, y, 0), then (x, y, 1). For example, the step set above is {(-1, 0, 0), (0, 1, 0), (1, 0, 0), (0, 0, 1), (0, 1, 1)}. The given numbers are the first coefficients in the expansion of G(t; 1, 1, 1).)

Discarding these cases from consideration, we are left with 35 different sequences whose generating series appear to be D-finite; among those, three appear to be algebraic. Their step sets are given in the appendix. We were not able to find an equation for the step set

[??] 1, 1, 4, 7, 28, 70, 280, 787, 3148, 9526, 38104, ... (A149080)

which is symmetric about all three axes, not even with 800 terms instead of 400. Also the step set

[??] 1, 1, 4, 13, 40, 136, 496, 1753, 6256, 22912, 85216, ... (A149424)

which enjoys a rotational symmetry about the middle line of the first octant, and which may be viewed as a three dimensional analogue of Kreweras's step set, appears to be non-D-finite, even when 800 terms are taken into account.

For walks in the quarter plane, it is conjectured in [29, Section 3] that D-finiteness is preserved under reversing arrows, i.e., the generating function for a step set [??] is D-finite if and only if the generating function for the step set [??]' is, when [??]' is obtained from [??] by reversing all arrows. Our computations do not suggest that this criterion also applies in 3D. Among the 134 sequences we found D-finite, there are 42 which correspond to step sets in [??] for whose counterpart in [??]' we were not able to find an equation. Among those, there are some which satisfy only very large equations, so that chances are that they remain D-finite upon reversing arrows, but with equations which are too large for us to find. Others satisfy quite small equations, for example the sequence A026378 whose step set is given above.

4.2 Algebraic Observations

As in the 2D case, it turns out that most of the minimal polynomials of the algebraic series define curves of genus 0, which therefore can be rationally parameterized. There are twelve cases of genus 1, these are elliptic curves. Some of them turn out to be isomorphic (over [bar.Q]). For example, those corresponding to the cases

[??] 1, 1, 4, 11, 32, 110, 360, 1163,4112, 14066, 47848, ... (A149232)

[??] 1, 2, 7, 27, 105, 426, 1787, 7590, 32633, 142152, 624659, ... (A150591)

[??] 1, 2, 10, 40, 176, 808, 3720, 17152, 81440, 384448, ... (A151023)

all have 1728 as j-invariant. They originate from the 2D Kreweras walks. Most interestingly, there are also three step sets originating from the 2D reverse Kreweras walks ([??]) for which the genus is 5 (!).

For the transcendental series, we could observe the same phenomenon as in 2D: all the operators factor as a product of a single irreducible operator of order two and several operators of order one. We therefore expect again that all these series admit a representation as a hypergeometric pullback. As an example, the generating function G(t; 1, 1, 1) of the sequence

[??] 1, 1, 2, 4, 10, 25, 70, 196, 588, 1764, ... (A005817)

can be written in the form

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

This representation was found by Mark van Hoeij. It is beyond the scope of the standard tools of Maple or Mathematica.

4.3 Analytic Observations

Also concerning asymptotics, similar remarks apply as in 2D. All coefficient sequences grow like [kappa][n.sup.[alpha]][[rho].sup.n] for some constants [kappa], [alpha], [rho], where [rho] is an integer or an algebraic number of degree 2 and [alpha] is a non-positive number. We have not gone through the laborious task of determining the constants [kappa].

Appendix

[TABLE 1 OMITTED]

[TABLE 2 OMITTED]

[TABLE 3 OMITTED]

Acknowledgments. We are grateful to the anonymous referees for their pertinent suggestions which helped to improve the presentation of this paper. We thank Pierre Nicodeme and Bruno Salvy, who carefully read a preliminary version and made several useful comments. We also thank Mark van Hoeij for computing the example pullback in 3D given in Section 4.2.

References

[1] Y. Andre. G-functions and geometry. Aspects of Mathematics, E13. Friedr. Vieweg & Sohn, Braunschweig, 1989.

[2] D. H. Bailey, J. M. Borwein, N. J. Calkin, R. Girgensohn, D. R. Luke, and V. H. Moll. Integer relation detection. In Experimental Mathematics in Action, pages 29-31. AK Peters, 2007.

[3] C. Banderier and P. Flajolet. Basic analytic combinatorics of directed lattice paths. Theoret. Comput. Sci., 281(1-2):37-80, 2002. Selected papers in honour of Maurice Nivat.

[4] B. Beckermann and G. Labahn. A uniform approach for the fast computation of matrix-type Pade approximants. SIAM J. Matrix Anal. Appl., 15(3):804-823, 1994.

[5] A. Bostan, F. Chyzak, G. Lecerf, B. Salvy, and E. Schost. Differential equations for algebraic functions. In C. W. Brown, editor, ISSAC'07: Proceedings of the 2007 international symposium on Symbolic and algebraic computation, pages 25-32. ACM Press, 2007.

[6] A. Bostan and M. Kauers. The complete generating function for Gessel walks is algebraic. In preparation, 2009.

[7] M. Bousquet-Melou. Counting walks in the quarter plane. In Mathematics and computer science, II (Versailles, 2002), Trends Math., pages 49-67. Birkhauser, Basel, 2002.

[8] M. Bousquet-Melou. Walks in the quarter plane: Kreweras' algebraic model. Ann. Appl. Probab., 15(2):1451-1491,2005.

[9] M. Bousquet-Melou and M. Mishna. Walks with small steps in the quarter plane. Available at http: //arxiv.org/abs/0810.4387,2008.

[10] M. Bousquet-Melou and M. Petkovsek. Linear recurrences with constant coefficients: the multivariate case. Discrete Math., 225(1-3):51-75, 2000. FPSAC'98 (Toronto, ON).

[11] R. Brak and A. J. Guttmann. Algebraic approximants: a new method of series analysis. J. Phys. A, 23(24):L1331-L1337, 1990.

[12] C. Brezinski. Algorithmes d'acceleration de la convergence. Editions Technip, Paris, 1978. Etude numerique, Collection Langages et Algorithmes de l'Informatique.

[13] A. Chambert-Loir. Theoremes d'algebricite en geometrie diophantienne (d'apres J.-B. Bost, Y. Andre, D. & G. Chudnovsky). Seminaire Bourbaki, 282(886):175-209, 2002.

[14] D. V. Chudnovsky and G. V. Chudnovsky. Applications of Pade approximations to Diophantine inequalities in values of G-functions. In Number theory (New York, 1983-84), volume 1135 of Lecture Notes in Math., pages 9-51. Springer, Berlin, 1985.

[15] M. Dettweiler and S. Reiter. On globally nilpotent differential equations. ArXiv:math/0605383, 2006.

[16] B. Dwork. Differential operators with nilpotent p-curvature. Amer. J. Math., 112(5):749-786, 1990.

[17] B. Dwork, G. Gerotto, and F. J. Sullivan. An introduction to G-functions, volume 133 of Annals of Mathematics Studies. Princeton University Press, Princeton, NJ, 1994.

[18] P. Flajolet. Analytic models and ambiguity of context-free languages. Theoretical Computer Science, 49(2-3):283-309, 1987.

[19] P. Flajolet and R. Sedgewick. Analytic Combinatorics. Cambridge University Press, 2009.

[20] S. Garoufalidis. G-functions and multisum versus holonomic sequences. Advances in Mathematics, 220(6):1945-1955, 2009.

[21] T. Honda. Algebraic differential equations. In Symposia Mathematica, Vol. XXIV (Sympos., INDAM, Rome, 1979), pages 169-204. Academic Press, London, 1981.

[22] N. M. Katz. Nilpotent connections and the monodromy theorem: Applications of a result of Turrittin. Inst. Hautes Etudes Sci. Publ. Math., 39:175-232, 1970.

[23] M. Kauers. Guessing handbook. Technical Report 09-07, RISC-Linz, 2009.

[24] M. Kauers, C. Koutschan, and D. Zeilberger. Proof of Ira Gessel's lattice path conjecture. Technical Report 2008-08, SFB F13, 2008.

[25] M. Kauers and D. Zeilberger. The quasi-holonomic ansatz and restricted lattice walks. J. Difference Equ. Appl., 14(10-11):1119-1126, 2008.

[26] D. Krammer. An example of an arithmetic Fuchsian group. J. Reine Angew. Math., 473:69-85, 1996.

[27] G. Kreweras. Sur une classe de problemes de denombrement lies au treillis des partitions des entiers. Cahiers du B.U.R.O., 6:5-105, 1965.

[28] C. Mallinger. Algorithmic manipulations and transformations of univariate holonomic functions and sequences. Master's thesis, J. Kepler University, Linz, August 1996.

[29] M. Mishna. Classifying lattice walks restricted to the quarter plane. J. Combin. Theory Ser. A, 116(2):460-177, 2009.

[30] S. Plouffe. Plouffe's inverter. http://pi.lacim.uqam.ca/.

[31] B. Salvy and P. Zimmermann. Gfun: a Maple package for the manipulation of generating and holonomic functions in one variable. ACM Transac. Math. Software, 20(2):163-177, 1994.

[32] H. Schmitt. Operators with nilpotent p-curvature. Proc. Amer. Math. Soc., 119(3):701-710, 1993.

[33] J. Wimp and D. Zeilberger. Resurrecting the asymptotics of linear recurrences. Journal of Mathematical Analysis and Applications, 111:162-176, 1985.

Alin Bostan (1) ([dagger]) and Manuel Kauers (2) ([double dagger])

(1) Algorithms Project, INRIA Paris-Rocquencourt, 78153 Le Chesnay, France

(2) Research Institute for Symbolic Computation, J. Kepler University Linz, Austria

([dagger]) Partially supported by the French National Agency for Research (ANR Project "Gecko") and the Microsoft Research-INRIA Joint Centre.

([double dagger]) Partially supported by the Austrian Science Foundation (FWF) grants P19462-N18 and P20162-N18.

(i) The usual definition is more general, the coefficients of S can be taken in an arbitrary algebraic number field. For our purposes it is sufficient and convenient to restrict to rational coefficients.

(ii) We have carried out our computations on various different machines whose main memory ranges from 8 Gb to 32 Gb and which are equipped with (multiple) processors all running at about 3GHz.

(iii) Modulo 5, the curvature matrix [M.sub.5](t) has [T.sup.2] as minimal polynomial.