Alin Bostan
email: alin.bostan@inria.fr

Mme Frédérique Bassino

M Olivier Bodini

Mme Mireille

M Florent Hivert

M Gilles Villard

COMPUTER ALGEBRA FOR LATTICE PATH COMBINATORICS

Keywords: Enumerative combinatorics, random walks in cones, lattice paths in the quarter plane, Gessel walks, generating functions, computer algebra, automated guessing, creative telescoping, diagonals, binomial sums, algebraic functions, D-finite functions, hypergeometric functions, elliptic integrals AMS subject classifications. Primary 05A10, 05A15, 05A16, 97N70, 33F10, 68W30, 14Q20; Secondary 33C05, 97N80, 13P15, 33C75, 12Y05, 13P05, 14Q20

Classifying lattice walks in restricted lattices is an important problem in enumerative combinatorics. Recently, computer algebra has been used to explore and to solve a number of difficult questions related to lattice walks. We give an overview of recent results on structural properties and explicit formulas for generating functions of walks in the quarter plane, with an emphasis on the algorithmic methodology.

A tandem-walk is a path in Z 2 taking steps from {↑, ←, } only. Show that, for any integer n ≥ 0, the following quantities are equal: (i) the number a n of tandem-walks of length n (i.e., using n steps), confined to the upper half-plane Z × N, that start and end at (0, 0); (ii) the number b n of tandem-walks of length n confined to the quarter plane N 2 , that start at (0, 0) and finish on the diagonal x = y.

For instance, for n = 3, this common value is a 3 = b 3 = 3, as shown below.

(i) (ii)

The problem establishes a rather surprising connection between tandem-walks in the lattice plane, submitted to two different kinds of constraints: the evolution domain of the walk, and its ending point. The domain constraint is weaker for the first family of walks, while the ending constraint is relaxed for the second family.

It appears that this problem is far from being trivial. Several solutions exist, but none of them is elementary. One of the main aims of the present text is to convince the reader that this problem (and many others with a similar flavor) can be solved with the help of a computer. More precisely, Computer Algebra tools, extensively described in the following sections, can be used to discover and to prove the following equalities [START_REF] Abel | OEuvres complètes. Tome II[END_REF] a 3n = b 3n = (3n)! n! 2 • (n + 1)! , and a m = b m = 0 if 3 does not divide m.

It goes without saying that such a simple and beautiful expression cannot be an element of chance. As it will turn out, closed forms are quite rare for this kind of enumeration problems. Nevertheless, even in absence of nice formulas, the structural properties of the corresponding enumeration sequences reflect the symmetries of the step set and of the evolution domain. Equation [START_REF] Abel | OEuvres complètes. Tome II[END_REF] shows that the sequences (a n ) and (b n ) are P-recursive, that is, they satisfy a linear recurrence with polynomial coefficients (in the index n). One of the messages that will emerge from the text is that this important property of the enumeration sequences is intimately related to the finiteness of a certain group, naturally attached to the step set {↑, ←, }.

1.2. General context: lattice paths confined to cones. Let us put the previous problem into a more general framework. Let d ≥ 1 be an integer (dimension), let S be a finite subset (called step set, or model) of vectors in Z d , and p 0 ∈ Z d (starting point). A S-path (or S-walk) of length n starting at p 0 is a sequence (p 0 , p 1 , . . . , p n ) of elements in the lattice Z d such that p i+1p i ∈ S for all 0 ≤ i < n. Let C be a cone of R d , that is a subset of R d such that r • v ∈ C for any v ∈ C and r 0, assumed to contain p 0 . We will be interested in the (exact and asymptotic) enumeration of S-walks confined to the cone C, and potentially subject to additional constraints.

Example 1. Consider the model S = {(1, 0), (-1, 0), (1, -1), (-1, 1)} (called the Gouyou-Beauchamps model) in dimension d = 2, with starting point p 0 = (0, 0) and with cone C = R 2 + (the quarter plane). The picture below displays the step set of the model (on the left), and a S-walk of length n = 17 confined to C (on the right). 00000000000 00000000000 00000000000 00000000000 00000000000 00000000000 00000000000 00000000000 00000000000 00000000000 00000000000 00000000000 11111111111 11111111111 11111111111 11111111111 11111111111 11111111111 11111111111 11111111111 11111111111 11111111111 11111111111 11111111111 (i, j) = (5, 1)

The main typical questions in this context are then the following:

• What is the number a n of n-step S-walks contained in C and starting at p 0 ?

• For fixed i ∈ C, what is the number a n;i of such walks that end at i?

• What is the nature of their generating functions A(t) = ∑ n a n t n and A(t; x) = ∑ n,i a n;i t n x i ?

As expected from the introductory example of tandem-walks, the answers to these questions are not simple, and heavily depend on the various parameters. The aim of this text is to provide a survey of recent results -notably classification results and closed form expressions-obtained using Computer Algebra.

Why count walks in cones?.

Lattice paths are fundamental objects in combinatorics. They have been studied at least since the second half of the 19th century, in connection with the ballot problem (see §1.4). Even earlier, embryonic occurrences (around 1650) are in Pascal's and Huygens' solutions of the so-called problem of division of the stakes (or, problem of points), and of the gambler's ruin problem, which motivated the beginnings of modern probability theory [START_REF] Edwards | Pascal's arithmetical triangle[END_REF][START_REF] Hald | A history of probability and statistics and their applications before 1750[END_REF][START_REF] Devlin | The unfinished game[END_REF]. Despite these historically important examples, the enumeration of lattice walks has long remained part of what may be called recreational mathematics. It is only in the late 1960s that their study really became an independent field of research, at the crossroads of pure and applied mathematics. Since then, various approaches have been progressively involved, separately or in interaction, in the study of lattice walks. These methods arise from various fields of classical mathematics (algebra, combinatorics, complex analysis, probability theory), and more recently from computer science. There are several reasons for the ubiquity of lattice walks, but the most solid one is that they encode several important classes of mathematical objects, in discrete mathematics (permutations, trees, words, urns, . . . ), in statistical physics (magnetism, polymers, . . . ), in probability theory (branching processes, games of chance, . . . ), in operations research (birth-death processes, queueing theory, . . . ). Therefore, many questions from all these various fields can be reduced to solving lattice path problems. For more motivations, the reader is referred to the introduction of [START_REF] Banderier | Basic analytic combinatorics of directed lattice paths[END_REF]. Nowadays, several books are entirely devoted to lattice paths and their applications [START_REF] Spitzer | Principles of random walk, The University Series in Higher Mathematics[END_REF][START_REF] Mohanty | Lattice path counting and applications[END_REF][START_REF] Narayana | Lattice path combinatorics with statistical applications[END_REF][START_REF] Doyle | Random walks and electric networks[END_REF][START_REF] Cohen | Analysis of random walks, vol. 2 of Studies in Probability, Optimization and Statistics[END_REF][START_REF] Weiss | Aspects and applications of the random walk, Random Materials and Processes[END_REF][START_REF] Fayolle | Random walks in the quarter plane[END_REF][START_REF] Woess | Random walks on infinite graphs and groups[END_REF][START_REF] Borovkov | Encyclopedia of Mathematics and its Applications[END_REF][START_REF] Lawler | Random walk: a modern introduction[END_REF][START_REF] Blanchard | Random walks and diffusions on graphs and databases[END_REF], and an international conference titled Lattice path combinatorics and applications is entirely devoted to this field. We recommend Humphreys' article [START_REF] Humphreys | A history and a survey of lattice path enumeration[END_REF] for a brief review of the history of lattice path enumeration and for a survey of the recent evolution of the field. Also, Krattenthaler's recent survey [START_REF] Krattenthaler | Lattice path enumeration[END_REF] is an excellent overview of various results and methods in lattice path enumeration.

The ballot problem and the reflection principle.

As mentioned before, the enumeration of lattice walks is an old topic. We want to illustrate this using Bertrand's ballot problem [START_REF] Bertrand | Solution d'un problème[END_REF][START_REF] André | Solution directe du problème résolu par M. Bertrand[END_REF]. The aim is not only to provide the flavor of a nice piece of combinatorial reasoning, but especially to introduce the so-called reflection principle, seemingly invented by Aebly and Mirimanoff [START_REF] Aebly | Démonstration du problème du scrutin par des considérations géométriques[END_REF][START_REF] Mirimanoff | À propos de l'interprétation géométrique du problème du scrutin[END_REF], which contains the roots of a systematic method for lattice walks, to be presented later, and based on the notion of group of a walk, see §1.18. Bertrand's problem is the following: Suppose that two candidates A and B are running in an election. If a votes are cast for A and b votes are cast for B, where a > b, then what is the probability that A stays (strictly) ahead of B throughout the counting of the ballots?

The problem admits an obvious lattice path reformulation. Let us call a Dyck path a walk in the lattice plane Z 2 , with step set S = {(1, 1), (1, -1)} = { , }, that starts at the origin. Then, the problem asks for the number of Dyck paths consisting of a upsteps and b downsteps such that no step ends on the x-axis. Let us call these good paths. Clearly, any such good path starts with a step from (0, 0) to [START_REF] Abel | OEuvres complètes. Tome II[END_REF][START_REF] Abel | OEuvres complètes. Tome II[END_REF], and finishes at the point T(a + b, ab). Instead of counting good paths, it is actually easier to count bad paths: these are Dyck paths consisting of a upsteps and b downsteps that touch the x-axis at least once. Now enters the crucial observation, based on a reflection argument (see the picture).

To any bad path one may bijectively attach an unconstrained path in Z 2 from (1, -1) to T by simply reflecting, with respect to the horizontal axis, the first portion of the walk, which lies strictly above the horizontal axis before touching it for the first time. Therefore, the number of good paths is exactly the difference between the unconstrained Dyck paths in Z 2 from (1, 1) to T(a + b, ab) and the unconstrained Dyck paths in Z 2 from (1, -1) to T(a + b, ab). Since unconstrained Dyck paths are simply counted by binomials, that number is:

a + b -1 a -1 - a + b -1 b -1 = a -b a + b a + b a ,
from which one directly deduces the answer (ab)/(a + b) to Bertrand's problem.

Observe that, when a = n + 1 and b = n, the number of good paths is the famous Catalan number

C n = 1 2n + 1 2n + 1 n + 1 = 1 n + 1 2n n ,
that counts a plethora of different combinatorial objects [START_REF] Catalan | Note sur une équation aux différences finies[END_REF][START_REF] Catalan | Solution nouvelle de cette question : Un polygône étant donné, de combien de manières peut-on le partager en triangles au moyen de diagonales ?[END_REF][START_REF] Stanley | Catalan numbers[END_REF].

There exists a second (non-strict) version of the problem, in which A has at least as many votes as B all along the counting. The reflection principle still applies, and the answer is 1b/(a + 1). More information, and historical background, on the ballot problem is provided in the articles [START_REF] Barton | Some aspects of the random sequence[END_REF][START_REF] Renault | Lost (and found) in translation: André's actual method and its application to the generalized ballot problem[END_REF].

Last, but not least, let us mention that a higher dimensional version of the reflection principle [START_REF] Grossman | Fun with lattice points[END_REF][START_REF] Zeilberger | André's reflection proof generalized to the many-candidate ballot problem[END_REF] can be used to solve the following generalization of the ballot problem: Assume there are d candidates in an election, say A 1 , . . . , A d , with each A i receiving a i votes. What is the probability that, throughout the counting of the ballots, A i has at least as many votes as A i+1 for all 1 ≤ i ≤ d -1? This amounts to counting paths in Z d from the origin to (a 1 , . . . , a d ) that use only unit positive steps (in the direction of some coordinate axis) and that are confined to the edge cone {x 1 ≥ x 2 ≥ • • • ≥ x d ≥ 0}. The natural setting for the most general version of the reflection principle is the one of reflection groups: it applies when the set of steps is left invariant by a Weyl group and the walks are confined to a corresponding Weyl chamber see [START_REF] Gessel | Random walk in a Weyl chamber[END_REF][START_REF] Grabiner | Random walks in Weyl chambers and the decomposition of tensor powers[END_REF] and [269, §10.18]. 1.5. Pólya's "promenade au hasard" / "Irrfahrt". Another old and famous result on lattice paths is Pólya's theorem [328, 329] * about the so-called drunkard walk in the d-dimensional integer lattice Z d . By definition, such a walk is a random path in Z d for the so-called simple model, or Pólya's model. After a busy night at the bar (some vertex of Z d ), a drunkard wishes to get home (another vertex of Z d ). Given his mental and physical state, he cannot do better than executing a random walk starting from the bar: at each tick of the clock he moves to one of the 2d neighbors of the current vertex, chosen uniformly at random. What is the probability that he ever reaches his destination? The interesting fact is that the long-term behavior of the drunkard's walk depends on the dimension d.

Theorem 2 [START_REF] Pólya | Über eine Aufgabe der Wahrscheinlichkeitsrechnung betreffend die Irrfahrt im Straßennetz[END_REF]. Consider the simple random walk on Z d . If d ∈ {1, 2}, then the walk returns to its starting position with probability 1 (the simple walk is recurrent). If d ≥ 3, then with positive probability, the walk never returns to its starting position (the simple walk is transient). † Several proofs exist for this classical result. Probably the most direct one [START_REF] Feller | An introduction to probability theory and its applications[END_REF]§XIV.7] is based on the observation that the probability for the d-dimensional drunkard to be back at the origin after 2n steps is equal to the (d -1)-folded sum

u (d) 2n = ∑ i 1 +•••+i d =n (2n)! (i 1 ! • • • i d !) 2 1 2d 2n .
Then some algebraic manipulations and Stirling's formula imply the asymptotic estimate u

(d) 2n = Θ(n -d/2
). On the other hand, it is not hard to see that the walk is transient if and only if the series ∑ n≥0 u (d) 2n converges, namely to a value m d which is the expected number of returns at the origin.

As a consequence of Theorem 2, if the drunkard lives in a 2-dimensional city, then he will eventually get home, even though possibly after a very long amount of time. But if, by misfortune, he lives in a 3-dimensional city, then the probability p 3 of return home will be less than 1. Pólya did not find a value for p 3 ; this was done later by McCrea and Whipple [START_REF] Mccrea | Random paths in two and three dimensions[END_REF] who showed that p 3 ≈ 0.34053. A beautiful exact formula for p 3 was found by Glasser and Zucker [START_REF] Glasser | Extended Watson integrals for the cubic lattices[END_REF], in terms of Euler's gamma function Γ(x) = ∞ 0 e -t t x-1 dt. It reads p 3 = 1 -1/m 3 , where

m 3 = √ 6 32π 3 Γ 1 24 Γ 5 24 Γ 7 24 Γ 11 24
≈ 1.516386060, see also [160, §2.3.5] and [START_REF] Borwein | Fast evaluation of the gamma function for small rational fractions using complete elliptic integrals of the first kind[END_REF][START_REF] Guttmann | Lattice Green's functions in all dimensions[END_REF][START_REF] Zucker | 70+ years of the Watson integrals[END_REF][START_REF] Koutschan | Lattice Green functions of the higher-dimensional face-centered cubic lattices[END_REF]. No similar closed-form expression is known for d ≥ 4, although it was proved [START_REF] Montroll | Random walks in multidimensional spaces, especially on periodic lattices[END_REF] that the probability of return p d equals 1 -1/m d , with

m d = d (2π) d π -π • • • π -π dx 1 • • • dx d d -cos x 1 -. . . -cos x d = ∞ 0 (I 0 (t/d)) d e -t dt,
where I 0 (t) is the modified Bessel function of the first kind I 0 (t) = ∑ k≥0 (t 2 /4) k k! 2 . A question closely related with Pólya's theorem will be discussed in §1.7.

Blending Experimental Mathematics and Computer Algebra in the service of lattice paths combinatorics.

The examples in §1.4 and §1.5 show that the study of lattice walks is an old field of research. The following sections will demonstrate that their exact and asymptotic enumeration is still a topical issue, with a lot of recent activity, new and exciting results, and many open questions. For instance, even when only restricting to articles published since 2000, and when only focusing to the case of walks confined to the quarter plane, one realizes that this particular case has received special attention, and much progress has been done by many recent contributors [START_REF] Choi | Counting lattice paths in restricted planes[END_REF][START_REF] Sulanke | Counting lattice paths by Narayana polynomials[END_REF][START_REF] Banderier | Limit laws for basic parameters of lattice paths with unbounded jumps[END_REF][START_REF] Banderier | Basic analytic combinatorics of directed lattice paths[END_REF][START_REF] Bousquet-Mélou | Counting walks in the quarter plane[END_REF][START_REF] Bousquet-Mélou | Walks in the quarter plane: a functional equation approach[END_REF][START_REF] Humphreys | Counting lattice paths with an infinite step set and special access[END_REF][START_REF] Bousquet-Mélou | Walks confined in a quadrant are not always D-finite[END_REF][START_REF] Humphreys | Counting lattice paths taking steps in infinitely many directions under special access restrictions[END_REF][START_REF] Bousquet-Mélou | Walks in the quarter plane: Kreweras' algebraic model[END_REF][START_REF] Niederhausen | Random walks in octants, and related structures[END_REF][START_REF] Bousquet-Mélou | Polynomial equations with one catalytic variable, algebraic series and map enumeration[END_REF][START_REF] Bernardi | Bijective counting of Kreweras walks and loopless triangulations[END_REF][START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF][START_REF] Kauers | The quasi-holonomic ansatz and restricted lattice walks[END_REF][START_REF] Ayyer | Towards a human proof of Gessel's conjecture[END_REF][START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF][START_REF] Kauers | Proof of Ira Gessel's lattice path conjecture[END_REF][START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF][START_REF] Mishna | Two non-holonomic lattice walks in the quarter plane[END_REF][START_REF] Böhm | Lattice path counting and the theory of queues[END_REF][START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF][START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF][START_REF] Fayolle | On the holonomy or algebraicity of generating functions counting lattice walks in the quarter-plane, Markov Process[END_REF][START_REF] Fayolle | Random walks in the quarter plane with zero drift: an explicit criterion for the finiteness of the associated group, Markov Process[END_REF][START_REF] Guillemin | Rare event asymptotics for a random walk in the quarter plane[END_REF][START_REF] Kurkova | Explicit expression for the generating function counting Gessel's walks[END_REF][START_REF] Kurkova | Random walks in (Z + ) 2 with non-zero drift absorbed at the axes[END_REF][START_REF] Fayolle | Some exact asymptotics in the counting of walks in the quarter plane[END_REF][START_REF] Kurkova | On the functions counting walks with small steps in the quarter plane[END_REF][START_REF] Raschel | Counting walks in a quadrant: a unified approach via boundary value problems[END_REF][START_REF] Sun | Proof of two conjectures of Petkovšek and Wilf on Gessel walks[END_REF][START_REF] Kung | Catalan lattice paths with rook, bishop and spider steps[END_REF][START_REF] Melczer | Singularity analysis via the iterated kernel method[END_REF][START_REF] Raschel | Walks in the quarter plane: analytic approach and applications[END_REF][START_REF] Raschel | Random walks in the quarter plane, discrete harmonic functions and conformal mappings[END_REF][START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF][START_REF] Duraj | Random walks in cones: the case of nonzero drift[END_REF][START_REF] Melczer | Enumerating lattice paths with symmetries through multivariate diagonals[END_REF][START_REF] Melczer | Asymptotic lattice path enumeration using diagonals[END_REF][START_REF] Albert | Permutations sortable by two stacks in parallel and quarter plane walks[END_REF][START_REF] Bostan | Compter les excursions sur un échiquier[END_REF][START_REF] Denisov | Random walks in cones[END_REF][START_REF] Fayolle | About a possible analytic approach for walks in the quarter plane with arbitrary big jumps[END_REF][START_REF] Fayolle | Random walks in the quarter-plane: advances in explicit criterions for the finiteness of the associated group in the genus 1 case, Markov Process[END_REF][START_REF] Kauers | Walks in the quarter plane with multiple steps[END_REF][START_REF] Kurkova | New steps in walks with small steps in the quarter plane: series expressions for the generating functions[END_REF][START_REF] Bacher | Continued classification of 3D lattice models in the positive octant[END_REF][START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF][START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF][START_REF] Bousquet-Mélou | Square lattice walks avoiding a quadrant[END_REF][START_REF] Bousquet-Mélou | An elementary solution of Gessel's walks in the quadrant[END_REF][START_REF] D'arco | Combinatorics meets potential theory[END_REF][START_REF] Garbit | On the exit time from a cone for random walks with drift[END_REF][START_REF] Melczer | Asymptotic lattice path enumeration using diagonals[END_REF][START_REF] Melczer | Asymptotics of lattice walks via analytic combinatorics in several variables[END_REF][START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF][START_REF] Bostan | A human proof of Gessel's lattice path conjecture[END_REF][START_REF] Courtiel | Weighted lattice walks and universality classes[END_REF][START_REF] Du | Infinite orders and non-D-finite property of 3-dimensional lattice walks[END_REF][START_REF] Kauers | Lattice walks in the octant with infinite associated groups[END_REF][START_REF] Mishna | Algorithmic approaches for lattice path combinatorics[END_REF]. And this is certainly not an exhaustive list. † As Feller says [185, p. 360], the statement "all roads lead to Rome" is justified in two dimensions.

The dominating point of view in these works is to develop uniform approaches, rather than ad-hoc solutions to a specific question. My personal bias is twofold: combine an experimental mathematics approach, as promoted in the beautiful and inspiring books by Borwein and collaborators [START_REF] Borwein | Experimentation in mathematics[END_REF][START_REF] Bailey | Experimental mathematics in action[END_REF][START_REF] Borwein | Mathematics by experiment[END_REF][START_REF] Borwein | The computer as crucible[END_REF], with modern tools from the Computer Algebra arsenal as described in the recent reference textbooks [START_REF] Zur Gathen | Modern computer algebra[END_REF][START_REF] Bostan | Algorithmes efficaces en calcul formel[END_REF], in order to conjecture and prove enumerative and asymptotic results for lattice paths.

Over the last three decades a fundamental shift has been operated in the way mathematics is practiced. As a consequence of the continued advance of computing power and of the unceasing availability of modern computational software, one can nowadays really take advantage of computer-aided research in order to solve significant and difficult mathematical problems. Our goal in this memoir is to overview computational approaches to discovery of new results in lattice path combinatorics. We entirely share Borwein's viewpoint that mathematical discovery through experimentation and the use of increasingly intelligent software is going to play an essential role in other fields of mathematics. 1.7. Another example, from the SIAM 100-Digit Challenge [START_REF] Trefethen | A hundred-dollar, hundred-digit challenge[END_REF][START_REF] Bornemann | The SIAM 100-digit challenge[END_REF]. In a 2002 SIAM News article [START_REF] Trefethen | A hundred-dollar, hundred-digit challenge[END_REF], L. N. Trefethen, head of the Numerical Analysis Group at Oxford University, proposed a contest which consisted of ten challenging problems in numerical computing. Each problem was stated in at most three simple sentences and had a single real number as a solution. The objective was to compute each number to as many digits of precision as possible. Scoring for the contest would be simple: each correct digit of the answer, up to ten per problem, would earn a single point. Trefethen warned that the problems were hard and indicated that he would be impressed if anyone managed to score even 50 points. Problem 6 in his list was about lattice walks in the plane, and appears to be related to Pólya's problem. Problem 6 (Biasing for a Fair Return) A flea starts at (0, 0) on the infinite two-dimensional integer lattice and executes a biased random walk: At each step it hops north or south with probability 1/4, east with probability 1/4 + ε, and west with probability 1/4ε. The probability that the flea returns to (0, 0) sometime during its wanderings is 1/2. What is ε?

As demonstrated in the wonderful book [START_REF] Bornemann | The SIAM 100-digit challenge[END_REF]Chap. 6], and in §3.2.1, Computer Algebra is able to conjecture and to prove the following formula

p( ) = 1 - A 2 • 2 F 1 1 2 , 1 2 1 2 √ 1 -16 2 A -1 , with A = 1 + 8 2 + 1 -16 2 , where 2 F 1 1 2 , 1 2 1 t = ∑ n≥0 2n n 2 t 16 n .
From this exact expression, it is easy to get the first 100 digits of the result ≈ 0.0619139544739909428481752164732121769996387749983 6207606146725885993101029759615845907105645752087861 . . . and actually millions of digits, if needed, in not more than a couple of seconds.

1.8. Two basic cones: the full space and a (rational) half-space. Let us now turn back to the general problem as stated in §1.2, using notion introduced in there. The simplest possible cone is the full space C = R d . In that case, the situation is very simple: the full generating function has the most basic structure, it is rational.

Theorem 3. If S ⊂ Z d and C = R d , then a n = |S| n , i.e. A(t) = ∑ n 0 a n t n = 1 1 -|S| t .
More generally:

A(t; x) = ∑ n,i a n;i x i t n = 1 1 -t ∑ s∈S x s .
The next case by increasing order of difficulty is when the cone is a half-space. The full generating function is not rational anymore, but nevertheless it still has a very important property: it is algebraic. Theorem 4. If S ⊂ Z d and if C is a rational half-space, then A(t; x) is algebraic, given by an explicit system of polynomial equations.

This result is due to Bousquet-Mélou and Petkovšek, see [START_REF] Bousquet-Mélou | Linear recurrences with constant coefficients: the multivariate case[END_REF]Theorem 13] and [START_REF] Bousquet-Mélou | Walks confined in a quadrant are not always D-finite[END_REF]Proposition 2]. Roots of it are in [START_REF] Petkovšek | The irrational chess knight[END_REF][START_REF] Petkovšek | Partial recurrences with constant coefficients: the algebraic case[END_REF]. The important particular case of 2D "generalized Dyck paths" had been treated before, see [START_REF] Gessel | A factorization for formal Laurent series and lattice path enumeration[END_REF][START_REF] Labelle | Generalized Dyck paths[END_REF][START_REF] Labelle | Langages de Dyck généralisés[END_REF][START_REF] Duchon | On the enumeration and generation of generalized Dyck words[END_REF]. The most basic illustration is provided by the ballot problem ( §1.4), for which A(t; 1) = ∑ n≥0 C n t n = (1 -√ 1 -4t)/(2t), see Example 5 below. The main ingredient in the proof of [START_REF] Bousquet-Mélou | Linear recurrences with constant coefficients: the multivariate case[END_REF] of Theorem 4, called the kernel method (terminology coined in [START_REF] Banderier | Basic analytic combinatorics of directed lattice paths[END_REF]), seems to belong to the "mathematical folklore". One source of this method, identified by Banderier and Flajolet in [26, p. 55], is Knuth's book [261, §2.2.1], more precisely his solutions to Exercises 4 and 11, which use a "new method for solving the ballot problem". Knuth's trick may have been better known at that time in probability theory, as suggested by its use in a more involved context [START_REF] Malyshev | Wiener-Hopf equations in the quarter-plane, discrete groups and automorphic functions[END_REF][START_REF] Malyshev | An analytical method in the theory of two-dimensional positive random walks[END_REF][START_REF] Flatto | Two queues in parallel[END_REF][START_REF] Fayolle | Solutions of functional equations arising in the analysis of two-server queueing models[END_REF][START_REF] Fayolle | Two coupled processors: the reduction to a Riemann-Hilbert problem[END_REF]. Various examples of its use in combinatorics are presented by Prodinger in [START_REF] Prodinger | The kernel method: a collection of examples[END_REF]. More historical notes on the origins of the kernel method can be found in [26, §2.2] and in [27, §1]. It is my feeling that the origins of the method amount at least to Kingman's article [START_REF] Kingman | Two similar queues in parallel[END_REF] in queueing theory, a reference that seems to have been previously overlooked. A very nice and powerful generalization of the kernel method is presented in [START_REF] Bousquet-Mélou | Polynomial equations with one catalytic variable, algebraic series and map enumeration[END_REF].

Example 5. Let us illustrate the kernel method on the simplest example, in relation with the ballot problem introduced in §1.4. Set S = {(1, 1), (1, -1)} = { , } and denote by M n,k be the number of S-walks in N 2 of length n that start at (0, 0) and end at vertical altitude k. Let M(x, y) = ∑ n,k M n,k x n y k . We will show that:

(a) M obeys the functional equation (y

-x(1 + y 2 )) • M(x, y) = y -x • M(x, 0). (b) M is algebraic, namely M(x, y) = √ 1 -4x 2 + 2xy -1 2x(y -x(1 + y 2 )) .
The starting point is an obvious recurrence relation, together with initial conditions, that translate the enumerative problem.

(2)

M n+1,k = M n,k-1 + M n,k+1 , M 0,0 = 1, M -1,k = M n,-1 = 0 for k, n ≥ 0.
Multiplying the recurrence relation by x n+1 y k+1 , and summing over n, k ∈ N yields

y • M(x, y) -∑ k≥0 M 0,k y k M(0,y) = 1 = y 2 x • M(x, y) + x • M -∑ n≥0 M n,0 x n M(x,0)
, which rewrites as the so-called kernel equation

(3) (y -x(1 + y 2 )) • M(x, y) = y -x • M(x, 0).
Observe that simple manipulations like setting y = 0 in (3) lead to tautologies. The kernel method consists in the following simple observation: let y 0 ∈ Q[[x]] be the power series root of K = yx(1 + y 2 ), the coefficient of M(x, y) in Eq. ( 3):

y 0 = 1 - √ 1 -4x 2 2x = x + x 3 + 2x 5 + 5x 7 + 14x 9 + • • • ∈ Q[[x]].
(One recognizes the generating function of Catalan numbers y 0 = ∑ n≥0 C n x 2n+1 .) Then, plugging y = y 0 into the kernel equation (3) delivers M(x, 0) = y 0 (x)/x. This provides an alternative, algebraic, proof of the (non-strict version of the) ballot problem. Finally, plugging back this value into (3) proves (b):

M(x, y) = y -y 0 K(x, y) = √ 1 -4x 2 + 2xy -1 2x(y -x(1 + y 2 )
) .

We will encounter more sophisticated uses of the kernel method in §2 and §3.

1.9. Lattice walks with small steps in the quarter plane. The next case by increasing level of complexity is the one of a cone obtained as the intersection of two half-spaces. Up to modifying the step set by a linear transformation, one may assume that the cone is the basic orthant C = R d + . This reduction is illustrated in the picture below, where the simple (Pólya) walks in the 2-dimensional cone of opening π/4 are put in bijection with the Gouyou-Beauchamps walks in the quarter plane.

(i, j) = (5, 1)

The power series expansions of many special functions in combinatorics and physics, including algebraic functions, are D-finite: they satisfy linear differential equations with polynomial coefficients, see §1.11 for definitions and main properties. For example, 60 % of the handbook [START_REF] Abramowitz | Handbook of mathematical functions with formulas, graphs, and mathematical tables[END_REF] describe D-finite functions.

That generating functions for walks constrained to evolve in an orthant need not be algebraic, and not even D-finite, was first observed by Bousquet-Mélou and Petkovšek in [START_REF] Bousquet-Mélou | Walks confined in a quadrant are not always D-finite[END_REF]. Preliminary results in this direction had been obtained by the same authors in [START_REF] Petkovšek | The irrational chess knight[END_REF][START_REF] Petkovšek | Partial recurrences with constant coefficients: the algebraic case[END_REF][START_REF] Bousquet-Mélou | Linear recurrences with constant coefficients: the multivariate case[END_REF]. The first model of walks in the quarter plane for which the generating function was proved to be non-D-finite [103, §3] is the socalled knight walks model: these are walks confined to N 2 that start from p 0 = (1, 1) and take their steps in S = {(2, -1), (-1, 2)}. This surprising result was the starting point of a massive classification effort, initiated by Mishna [START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF][START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF], intensified in a germinal work by Bousquet-Mélou and Mishna [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF], and continued by many researchers [START_REF] Kauers | The quasi-holonomic ansatz and restricted lattice walks[END_REF][START_REF] Ayyer | Towards a human proof of Gessel's conjecture[END_REF][START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF][START_REF] Kauers | Proof of Ira Gessel's lattice path conjecture[END_REF][START_REF] Mishna | Two non-holonomic lattice walks in the quarter plane[END_REF][START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF][START_REF] Kurkova | On the functions counting walks with small steps in the quarter plane[END_REF][START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF][START_REF] Kurkova | New steps in walks with small steps in the quarter plane: series expressions for the generating functions[END_REF]. The rest of this section is devoted to tell the story of this classification, with a viewpoint towards computerized proofs.

Before restricting our attention to the special but important case of walks with small steps in the quarter plane, let us mention two general criteria that contain sufficient conditions for D-finiteness of the full generating function A(t; x). One was obtained by Bousquet-Mélou in [94, §3]. (A combinatorial proof for the particular case of the length generating function A(t; 1) was given in [103, §2].) Theorem 6. Let C = R 2 + and let S ⊂ Z × {-1, 0, 1} be symmetric with respect to the horizontal axis. Then A(t; x) is D-finite, given by an explicit system of linear differential equations.

The other criterion, whose precise statement is too involved to be given here, was already mentioned in §1.4 in connection with the reflection principle. Its underlying idea (an algebraic version of the reflection principle) was discovered independently by Gessel and Zeilberger [START_REF] Gessel | Random walk in a Weyl chamber[END_REF] and Biane [START_REF] Biane | Minuscule weights and random walks on lattices[END_REF]. Roughly, the result asserts the following: if the set of steps is left invariant by a finite Weyl group, if the cone where the walks are confined to is a corresponding Weyl chamber and if no allowed step can traverse the boundary of the cone, then the generating function A(t; x) is D-finite. The precise assumptions can be found in [START_REF] Gessel | Random walk in a Weyl chamber[END_REF] and in [START_REF] Krattenthaler | Lattice path enumeration[END_REF]Th. 10.18.3]. The criterion then follows by combining [START_REF] Gessel | Random walk in a Weyl chamber[END_REF]Th. 3] with results on D-finiteness of positive parts and constant terms such as [START_REF] Lipshitz | The diagonal of a D-finite power series is D-finite[END_REF] (see also §3 of this document).

From now on, we focus on small-step walks (or, nearest-neighbor walks) in the quarter plane. These are walks in the lattice Z 2 , confined to the cone C = R 2 + (we will often say confined to N 2 ), that start at p 0 = (0, 0) and use steps in a model S which is a fixed subset of { , ←, , ↑, , →, , ↓}.

An example of a small-step walk for the model S = { , ←, ↑, →, , ↓}, with length n = 45 and ending point (i, j) = (14, 2), is depicted below.

S =

Let us denote by f n;i,j the number of walks of length n ending at (i, j). The full counting sequence ( f n;i,j ) n,i,j admits several interesting specializations:

• f n;0,0 , the number of walks of length n returning to origin ("excursions");

• f n = ∑ i,j≥0 f n;i,j , the number of walks with prescribed length n. As customary in combinatorics, to these enumeration sequences one attaches (univariate, or multivariate) power series, namely the complete generating function

F S (t; x, y) = ∞ ∑ n=0 ∞ ∑ i,j=0 f n;i,j x i y j t n ∈ Q[x, y][[t]],
and its corresponding univariate specializations:

• F S (t; 0, 0), the generating function of excursions;

• F S (t; 1, 1) = ∑ n≥0 f n t n , the length generating function;

• F S (t; 1, 0), resp. F S (t; 0, 1), the generating function of walks ending on the horizontal, resp. vertical, axis, also called boundary returns; • "F S (t; 0, ∞)" := x 0 F S (t; x, 1/x), the generating function of walks ending on the diagonal x = y of N 2 , also called diagonal returns. The general questions addressed in §1.2 specialize to the quarter-plane setting as follows: Given the model S, what can be said about the generating function F S (t; x, y), resp. about the counting sequence ( f n;i,j ) i,j,n , and their specializations? More precise sub-questions concern structures, explicit forms and asymptotics:

• Structures: is F S algebraic? Is it D-finite? None of them?

• Explicit forms: do F S (t; x, y) and ( f n;i,j ) i,j,n admit closed-form expressions?

• Asymptotics: what is the behavior of ( f n;0,0 ) n , and ( f n ) n when n → ∞? The emphasis will be put on how Computer Algebra can be used to give computational answers to these questions.

1.10. Small-step models of interest. Among the 2 8 models S ⊆ {-1, 0, 1} 2 \ {(0, 0)}, some are trivial (e.g., if S ⊆ { , ←, , , ↓}, then F S (t; x, y) ≡ 1), others are intrinsic to the half-plane (therefore F S (t; x, y) is algebraic, cf. Theorem 4), others come in pairs by diagonal symmetry (if S and S are symmetric with respect to the diagonal of N 2 , then F S (t; x, y) ≡ F S (t; y, x)), see Fig. 1. After discarding these cases, Bousquet-Mélou and Mishna [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF] found that there are exactly 79 interesting distinct models of small-step walks in the quarter plane. They are represented in Fig. 2, and are grouped in two classes: 74 nonsingular models (or genus-1 models in the terminology of [START_REF] Fayolle | Random walks in the quarter plane[END_REF]) and 5 singular models (or genus-0 models). Singular models are the ones for which walks never return to the origin, that is for which the excursions generating function is trivial F(t; 0, 0) ≡ 1. Among the 79 models, there are "special" ones, that are considered interesting enough and were enough studied to deserve names: Pólya:

; Kreweras: ; Gessel:

; Gouyou-Beauchamps: ; King: ; Tandem: . One objective is then to understand and classify all these 79 models according to the structural properties of their generating functions.

1.11. Classification of power series. Before stating the main results, we still need a few definitions on (univariate and multivariate) power series.

Definition 7. Let S(t) = ∑ ∞ n=0 s n t n be a power series in Q[[t]]. Then, S(t) is called • algebraic if it is a root of a non-trivial polynomial P ∈ Q[t, T], i.e., P t, S(t) = 0; • transcendental if it is not algebraic; • D-finite (or holonomic) if it is satisfies a non-trivial linear differential equation p r (t)S (r) (t) + • • • + p 0 (t)S(t) = 0 with polynomial coefficients p i (t) ∈ Q[t]; • hypergeometric if its coefficients sequence (s n ) n satisfies a non-trivial linear homo- geneous recurrence of order 1 with polynomial coefficients in Q[n].
A very important class of hypergeometric series is that of Gauss hypergeometric functions 2 F 1 with parameters a, b, c ∈ Q, c / ∈ -N, defined by

2 F 1 a b c t = ∞ ∑ n=0 (a) n (b) n (c) n t n n! , where (x) n = x(x + 1) • • • (x + n -1)
is the Pochhammer symbol. This notion admits an obvious extension to the so-called generalized hypergeometric function p F q depending on p + 1 rational parameters appearing in the top Pochhammer symbols, and on q rational parameters on the bottom. For example,

3 F 2 a b c d e t = ∞ ∑ n=0 (a) n (b) n (c) n (d) n (e) n t n n!
, where a, b, c, d, e ∈ Q and d, e / ∈ -N.

The way these three important classes of power series (algebraic, D-finite, hypergeometric) are connected is illustrated in Fig. 3.

That hypergeometric series are D-finite is an immediate consequence of the simple fact that coefficient sequences of D-finite series are exactly P-recursive sequences, satisfying linear recurrences with polynomial coefficients [START_REF] Stanley | Differentiably finite power series[END_REF].

That algebraic series are D-finite has been observed in 1827 by Abel [1, p. 287]. Cockle [START_REF] Cockle | On transcendental and algebraic solution[END_REF] gave an algorithm for the computation of such a differential equation of the minimal possible order, that Harley [START_REF] Harley | On the theory of the transcendental solution of algebraic equations[END_REF] called differential resolvent. The method was then rediscovered by Tannery [372, §17], see also [212, §2.4]. One of the applications of these differential equations is the efficient power series expansions of algebraic series: a linear differential equation translates into a linear recurrence, with the consequence that the number of operations required to compute the first N coefficients grows only linearly with N. This method has been popularized in the combinatorics community by Comtet [START_REF] Comtet | Calcul pratique des coefficients de Taylor d'une fonction algébrique[END_REF] and studied from the complexity point of view by Chudnovsky and Chudnovsky [START_REF] Chudnovsky | On expansion of algebraic functions in power and Puiseux series[END_REF][START_REF] Chudnovsky | On expansion of algebraic functions in power and Puiseux series[END_REF], and more recently in [START_REF] Bostan | Differential equations for algebraic functions[END_REF].

Finally, understanding power series that are simultaneously algebraic and hypergeometric is an old and difficult question. Fuchs asked in 1866 [START_REF] Fuchs | Zur Theorie der linearen Differentialgleichungen mit veränderlichen Coefficienten[END_REF] for a classification of all Gauss hypergeometric functions 2 F 1 a b c t that are algebraic. Fuchs' question was solved in 1873 by Schwarz [START_REF] Schwarz | Über diejenigen Fälle, in welchen die Gaußische hypergeometrische Reihe einer algebraische Funktion ihres vierten Elementes darstellt[END_REF], who showed using geometric arguments (sphere tilings by spherical triangles) that, up to some normalization of the parameters, and apart from an explicitly given finite number of sporadic cases,

2 F 1 r 1 -r 1 2 t = cos((1 -2r) • arcsin( √ t)) √ 1 -t , r ∈ Q
is the only family of algebraic 2 F 1 functions. Building on work by Eisenstein [START_REF] Eisenstein | Über eine allgemeine Eigenschaft der Reihen-Entwicklungen aller algebraischen Funktionen[END_REF][START_REF] Heine | Der Eisensteinsche Satz über Reihen-Entwickelung algebraischer Functionen[END_REF], Landau [START_REF] Landau | Eine Anwendung des Eisensteinschen Satzes auf die Theorie der Gaussschen Differentialgleichung[END_REF][START_REF] Landau | Über einen zahlentheoretischen Satz und seine Anwendung auf die hypergeometrische Reihe[END_REF] and Stridsberg [START_REF] Stridsberg | Sur le théorème d'Eisenstein et l'équation différentielle de Gauss[END_REF], Errera [START_REF] Errera | Zahlentheoretische Lösung einer functionentheoretischen Frage[END_REF]

, . . . , b k-1 , b k = 1} be two subsets of Q, assumed disjoint modulo Z. Let D be their common denominator. Then k F k-1 a 1 a 2 • • • a k b 1 • • • b k-1 t
is algebraic if and only if {e 2iπra j , j ≤ k} and {e 2iπrb j , j < k} interlace on the unit circle for all 1 ≤ r < D with gcd(r, D) = 1.

For instance, the following hypergeometric function [START_REF] Rodriguez-Villegas | Integral ratios of factorials and algebraic hypergeometric functions[END_REF], arising from Chebychev's work on the distribution of primes numbers [START_REF] Tchebichef | Mémoire sur les nombres premiers[END_REF] ∑ n (30n)!n! (15n)!(10n)!(6n)! t n = 8 F 7 2 14 3 9 5 5 t is an algebraic power series. Indeed, for all 1 ≤ r < 30 with gcd(r, 30) = 1, one obtains the picture in Fig. 4, where red circles that correspond to upper parameters of the 8 F 7 , are interlaced with blue circles that correspond to lower parameters.

Similar definitions for algebraicity and D-finiteness apply to multivariate power series. For instance, S ∈ Q[[x, y, t]] is algebraic if it is the root of a non-trivial polynomial P ∈ Q[x, y, t, T], and it is D-finite if the set of all partial derivatives of S spans a finite-dimensional vector space over Q(x, y, t), in other words if S satisfies a system

i a i (t, x, y) ∂ i S ∂x i = 0, ∑ i b i (t, x, y) ∂ i S ∂y i = 0, ∑ i c i (t, x, y) ∂ i S ∂t i = 0.
As in the univariate case, multivariate algebraic series are D-finite [START_REF] Lipshitz | D-finite power series[END_REF].

The concept of hypergeometric series also admits extensions to several variables, but they are beyond the scope of the present text. One such generalization was introduced around 1988 by Gel'fand, Kapranov and Zelevinsky [START_REF] Gel | Holonomic systems of equations and series of hypergeometric type[END_REF]200,[START_REF] Gel | Hypergeometric functions and toric varieties[END_REF][START_REF] Gel | Generalized Euler integrals and A-hypergeometric functions[END_REF][START_REF] Dwork | Hypergeometric series[END_REF][START_REF] Stienstra | GKZ hypergeometric structures[END_REF] and is known as GKZ-hypergeometric functions, or A-hypergeometric functions. Let us just mention that Beukers [START_REF] Beukers | Algebraic A-hypergeometric functions[END_REF] obtained a characterization of the class of algebraic GKZ-hypergeometric functions, that extends the interlacing criterion from [START_REF] Beukers | Monodromy for the hypergeometric function n F n-1[END_REF].

1.12. Kreweras' walks. An interesting model in the world of quarter-plane walks is Kreweras' model S = {↓, ←, }. It is related to a version of the threecandidate ballot problem, more difficult than the one mentioned at the end of §1.4. Let A, B, C be candidates in an election, that receive a, b, c votes respectively. What is the probability p(a, b, c) that, throughout the counting of the ballots, A has at least as many votes as B and at least as many votes as C? This amounts to counting paths in Z 3 from the origin to (a, b, c) that use only unit positive steps and that are confined to the cone {x 1 ≥ max(x 2 , x 3 ) ≥ 0} of Z 3 . It appears that the reflection principle does not apply here, contrary to the case of the edge cone {x 1 ≥ x 2 ≥ x 3 ≥ 0}.

Equivalently, the question amounts to counting paths in the quarter plane for the model S = {↓, ←, }. In a long paper, Kreweras [START_REF] Kreweras | Sur une classe de problèmes de dénombrement liés au treillis des partitions des entiers[END_REF] obtained a closedformula for p(a, b, c) as a binomial double-sum:

p(a, b, c) = 1 - b + c a + 1 + 1 (a + 1)(a + 2) b ∑ i=1 c ∑ j=1 b i c j 2i + 2j -2 2i -1 i + j + a a + 2 ,
which simplifies to P(a, b, 0) = 1b/(a + 1) for the two-candidate ballot problem (cf. §1.4), and to a simple formula in the special case c = a:

(4) p(a, b, a) = 2 2b+1 a! (a -b)! 2 (2a -2b + 1)! (2a + 2)! .
The same problem was considered independently by Flatto and Hahn [START_REF] Flatto | Two parallel queues created by arrivals with two demands I[END_REF] in an applied probabilistic context (double queue that arises when arriving customers simultaneously place two demands handled independently by two servers). As a consequence of Eq. ( 4), Kreweras obtained the following result, which was reproved using various methods in [START_REF] Kreweras | Solution of an enumerative problem connected with lattice paths[END_REF][START_REF] Niederhausen | The ballot problem with three candidates[END_REF][START_REF] Gessel | A probabilistic method for lattice path enumeration[END_REF][START_REF] Bousquet-Mélou | Counting walks in the quarter plane[END_REF][START_REF] Bousquet-Mélou | Walks in the quarter plane: Kreweras' algebraic model[END_REF][START_REF] Bernardi | Bijective counting of Kreweras walks and loopless triangulations[END_REF][START_REF] Kauers | The quasi-holonomic ansatz and restricted lattice walks[END_REF][START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF]. The last two references in this list provide two different computer-aided proofs. In what follows, we denote by K(t; x, y) = F S (t; x, y) the full generating function for Kreweras walks S = {↓, ←, } in the quarter plane, and by K(t; 0, 0) the generating function for Kreweras excursions. Theorem 9 (Kreweras, [START_REF] Kreweras | Sur une classe de problèmes de dénombrement liés au treillis des partitions des entiers[END_REF]). The generating function K(t; 0, 0) is equal to

E T c E ' ' T E E c E E T d d s d d E ' E ' d d s d d E E E d d d d d d s T c E ' ' T E E T E T E T © E ' E ' E E
(5) 3 F 2 1/3 2/3 1 3/2 2 27 t 3 = ∞ ∑ n=0 4 n ( 3n n ) (n + 1)(2n + 1) t 3n = 1 + 2t 3 + 16t 6 + 192t 9 + • • • .
As a corollary of Theorem 9, the results in §1.11 (e.g., Theorem 8) imply that K(t; 0, 0) is an algebraic power series. In fact, much more is true: Theorem 10 ( [START_REF] Flatto | Two parallel queues created by arrivals with two demands I[END_REF][START_REF] Gessel | A probabilistic method for lattice path enumeration[END_REF][START_REF] Bousquet-Mélou | Walks in the quarter plane: Kreweras' algebraic model[END_REF]). The full generating function K(t; x, y) for the Kreweras walks is algebraic.

In §2 we will sketch a computer-aided proof of this result [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF] based on the guess-and-prove paradigm.

Gessel's walks.

Probably the most difficult model of walks in the quarter plane is Gessel's model S = { , , ←, →}. In 2001, Ira Gessel formulated, in private conversations with colleagues (including Mireille Bousquet-Mélou, Doron Zeilberger and Guoce Xin), two conjectures equivalent to the following statements: Conjecture 1. The generating function G(t; 0, 0) of Gessel excursions is equal to

3 F 2 5/6 1/2 1 5/3 2 16t 2 = ∞ ∑ n=0 (5/6) n (1/2) n (5/3) n (2) n (4t) 2n = 1 + 2t 2 + 11t 4 + 85t 6 + • • • . Conjecture 2.
The full generating function G(t; x, y) is not D-finite.

Here, as for the Kreweras walks, we denoted by G(t; x, y) = F S (t; x, y) for S = { , , ←, →} the full generating function for Gessel walks in the quarter plane, and by G(t; 0, 0) the generating function for Gessel excursions.

The genesis of Gessel's conjectures is related to his interest in finding examples of cones in Z 2 for which the generating functions for the simple (Pólya's) walk would admit nice formulas. As discussed in §1.5, Pólya [START_REF] Pólya | Über eine Aufgabe der Wahrscheinlichkeitsrechnung betreffend die Irrfahrt im Straßennetz[END_REF] first observed that there are exactly ( 2n n ) 2 simple excursions of length 2n in the plane Z 2 , and that the full generating function is rational in that case. Still for the Pólya model, but now restricted to the half plane, resp. to the quarter plane, Arquès [START_REF] Arquès | Dénombrements de chemins dans R 2 soumis à contraintes[END_REF] proved that excursions of length 2n are counted by nice formulas: ( 2n+1 n )C n for Z × N, and C n C n+1 for N 2 . Concerning the nature of the full generating function, it is algebraic for the cone Z × N [START_REF] Bousquet-Mélou | Linear recurrences with constant coefficients: the multivariate case[END_REF], and D-finite for the cone N 2 [START_REF] Bousquet-Mélou | Counting walks in the quarter plane[END_REF]. Gouyou-Beauchamps [START_REF] Gouyou-Beauchamps | Chemins sous-diagonaux et tableaux de Young[END_REF] found a similar formula C n C n+2 -C 2 n+1 for the number of simple excursions of length 2n in the cone with angle 45 • (the first octant). The generating function for this cone is again D-finite [START_REF] Gessel | Random walk in a Weyl chamber[END_REF]. It was thus natural to consider the cone with angle 135 • , and this is what Gessel did. See [START_REF] Bostan | Compter les excursions sur un échiquier[END_REF] for more historical details. Moreover, G(t; x, y) is completely characterized by the functional equation ( 6): it is its unique solution in Q[x, y][[t]], and even in the ring Q[[x, y, t]]. Therefore, the task is simply to solve equation [START_REF] Aeppli | À propos de l'interprétation géométrique du problème du scrutin[END_REF].

-t 1 x + 1 x 1 y G(t; 0, y) -t 1 xy (G(t; x, 0) -G(t; 0, 0)). ( 6
)
Similarly, to any of the 79 models introduced in §1.10 is attached a very similar functional equation. Again, this equation merely reflects a step-by-step construction of quarter-plane walks, and is based on the most elementary decomposition: a walk is either the empty walk, or it is a shorter walk followed by a permissible step. This observation is naturally translated into a generating function equation using the inventory χ S (x, y) := ∑ (i,j)∈S x i y j , and the kernel K S (t; x, y) = xy(1tχ S (x, y)). Note that for a non-trivial model with small steps the kernel is a polynomial. The decomposition is translated into the kernel equation (we omit the subscript S): (7) K(t; x, y)F(t; x, y) = xy + K(t; x, 0)F(t; x, 0) + K(t; 0, y)F(t; 0, y) -K(t; 0, 0)F(t; 0, 0).

Observe that the last term of the right-hand side occurs only if the step belongs to the model S.

Following Zeilberger's terminology [START_REF] Zeilberger | The umbral transfer-matrix method. I. Foundations[END_REF], the variables x and y are called catalytic for equation [START_REF] Albert | Permutations sortable by two stacks in parallel and quarter plane walks[END_REF]. (This means that one cannot simply set x = 0 or y = 0 in the equation to solve for F(t; x, 0) and F(t; 0, y) first.) The number of catalytic variables is related to the number of constraints imposed to the walk. The case of kernel equations with a single catalytic variable corresponds to uni-directional walks and it is well-understood, the solutions being always algebraic [START_REF] Bousquet-Mélou | Linear recurrences with constant coefficients: the multivariate case[END_REF], see Theorem 4.

Classifying lattice walks in the quarter plane thus amounts to solving 79 such equations. In the remaining part of Section 1 we describe several classes of results in this direction that have been obtained using Computer Algebra tools.

1.15. Main results (I): algebraicity of Gessel walks. After an almost successful attempt in [START_REF] Kauers | The quasi-holonomic ansatz and restricted lattice walks[END_REF], Gessel's first conjecture was finally solved in 2009 by Kauers, Koutschan and Zeilberger in [START_REF] Kauers | Proof of Ira Gessel's lattice path conjecture[END_REF] using an extension of the guess-and-prove approach described in [START_REF] Kauers | The quasi-holonomic ansatz and restricted lattice walks[END_REF]. Theorem 11 ([254]). G(t; 0, 0

) = 3 F 2 5/6 1/2 1 5/3 2 16t 2 .
This result implies in particular that G(t; 0, 0) is D-finite, but has no immediate implications concerning the D-finiteness of G(t; x, y). It came as a total surprise when Bostan and Kauers [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF] proved that Gessel's second conjecture was false.

Theorem 12 ([85]). The generating function G(t; x, y) for Gessel walks is algebraic.

Prior to this result, even the algebraicity of G(t; 0, 0) had been overlooked, even though the classical results recalled in §1.11 obviously apply. For instance, because of the alternative representation [START_REF] Alexanderson | The random walks of George Pólya, MAA Spectrum[END_REF] 3

F 2 5/6 1/2 1 5/3 2 16t 2 = 1 t 2 1 2 2 F 1 -1/6 -1/2 2/3 16t 2 - 1 2 ,
it is clear that algebraicity of G(t; 0, 0) could have been decided using Schwarz's classification, but it appears that, quite strangely, nobody recognized that the parameters (-1/6, -1/2; 2/3) actually fit to Case III of Schwarz's table [START_REF] Schwarz | Über diejenigen Fälle, in welchen die Gaußische hypergeometrische Reihe einer algebraische Funktion ihres vierten Elementes darstellt[END_REF].

The original discovery and proof of Theorem 12 was computer-driven, and used a guess-and-prove approach, based on Hermite-Padé approximants. This will be explained in more details in §2. Note that as a byproduct of this proof, an estimate on the size of the minimal polynomial of G(t; x, y) has been given: according to [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF], that minimal polynomial has more than 10 11 terms when written in dense (expanded) form, for a total size of ≈ 30 Gb (!) Several human proofs of Theorem 12 have been discovered since the publication of [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF]: the first one used complex analysis [START_REF] Bostan | A human proof of Gessel's lattice path conjecture[END_REF], the second one was purely algebraic [START_REF] Bousquet-Mélou | Square lattice walks avoiding a quadrant[END_REF], and the more recent one is probably the most elementary [START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF][START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF]. These proofs also contain a proof of Theorem 11.

Main results (II):

Explicit form for G(t; x, y). An interesting consequence of Theorem 12 is the following result, which contains a closed-formula for the full generating function G(t; x, y) of Gessel walks [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF]. Then G(t; x, y) is equal to

Theorem 13 ([85]). Let V = 1 + 4t 2 + 36t 4 + 396t 6 +• • • be the unique root in Q[[t]] of (V -1)(1 + 3/V) 3 = (16t) 2 , let U = 1 + 2t 2 + 16t 4 + 2xt 5 + 2(x 2 + 83)t 6 + • • • be the unique root in Q[x][[t]] of x(V -1)(V + 1)U 3 -2V(3x + 5xV -8Vt)U 2 -xV(V 2 -24V -9)U + 2V 2 (xV -9x -8Vt) = 0,
and let W = t 2 + (y + 8)t 4 + 2(y 2 + 8y + 41)t 6 + • • • be be the unique root in Q[y][[t]] of y(1 -V)W 3 + y(V + 3)W 2 -(V + 3)W + V -1 = 0. OEIS S Pol
64(U(V+1)-2V)V 3/2 x(U 2 -V(U 2 -8U+9-V)) 2 -y(W-1) 4 (1-Wy)V -3/2 t(y+1)(1-W)(W 2 y+1) 2 (1 + y + x 2 y + x 2 y 2 )t -xy - 1 tx(y + 1)
.

Again, the original discovery and proof of this result was computer-driven. During the computerized proof, a few other remarkable facts have been noticed, namely that G(t; x, y) can be expressed using nested radicals; for instance the length generating function G(t; 1, 1)

= 1 + 2t + 7t 2 + 21t 3 + 78t 4 + • • • reads G(t; 1, 1) = - 1 2t + √ 3 6t H(t) + 16t(2t + 3) + 2 (1 -4t) 2 H(t) -H(t) 2 + 3 ,
where

H(t) = 1 + 4t 1/3 (1 + 4t) 2/3 /(1 -4t) 4/3 .
Actually, the proof uses the minimal polynomials for G(t; x, 0) and G(t; 0, y) that were guessed and proved during the algebraicity proof. A striking feature of Theorem 13 is the relative simplicity of the closed-form expression, especially when compared to the size of the minimal polynomial of G(t; x, y). As in the case of Theorem 12, the result in Theorem 13 admits several recent human proofs [START_REF] Bostan | A human proof of Gessel's lattice path conjecture[END_REF][START_REF] Bousquet-Mélou | Square lattice walks avoiding a quadrant[END_REF][START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF][START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF].

Main results (III): Models with D-Finite length generating function.

The computer-driven approach that allowed Bostan and Kauers [START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF] to discover and prove the properties of the puzzling generating function for Gessel walks was used as soon as 2008 by the same authors to provide a (conjecturally) exhaustive list of models having (conjecturally) D-finite and algebraic generating functions. That resulted in an experimental classification synthesized in Fig. 7, which displays 23 models of walks in the quarter plane for which the length generating function F(t; 1, 1) was conjectured to be D-finite. The computerized discovery used again a guess-and-prove method, based on Hermite-Padé approximation. Details will be OEIS S algebraic? asymptotics OEIS S algebraic? asymptotics

1 A005566 N 4 π 4 n n 13 A151275 N 12 √ 30 π (2 √ 6) n n 2 2 A018224 N 2 π 4 n n 14 A151314 N √ 6λµC 5/2 5π (2C) n n 2 3 A151312 N √ 6 π 6 n n 15 A151255 N 24 √ 2 π (2 √ 2) n n 2 4 A151331 N 8 3π 8 n n 16 A151287 N 2 √ 2A 7/2 π (2A) n n 2 5 A151266 N 1 2 3 π 3 n n 1/2 17 A001006 Y 3 2 3 π 3 n n 3/2 6 A151307 N 1 2 5 2π 5 n n 1/2 18 A129400 Y 3 2 3 π 6 n n 3/2 7 A151291 N 4 3 √ π 4 n n 1/2 19 A005558 N 8 π 4 n n 2 8 A151326 N 2 √ 3π 6 n n 1/2 9 A151302 N 1 3 5 2π 5 n n 1/2 20 A151265 Y 2 √ 2 Γ(1/4) 3 n n 3/4 10 A151329 N 1 3 7 3π 7 n n 1/2 21 A151278 Y 3 √ 3 √ 2Γ(1/4) 3 n n 3/4 11 A151261 N 12 √ 3 π (2 √ 3) n n 2 22 A151323 Y √ 23 3/4 Γ(1/4) 6 n n 3/4 12 A151297 N √ 3B 7/2 2π (2B) n n 2 23 A060900 Y 4 √ 3 3Γ(1/3) 4 n n 2/3 A = 1 + √ 2, B = 1 + √ 3, C = 1 + √ 6, λ = 7 + 3 √ 6, µ = 4 √ 6-1 19
Figure 8. Models with D-Finite length generating function F S (t; 1, 1): asymptotics of f n = [t n ]F(t; 1, 1). For models 11, 13 and 15, estimates only hold for even n; for odd n, the constants change into 18 π , 144 √ 5 and 32 π [START_REF] Melczer | Asymptotics of lattice walks via analytic combinatorics in several variables[END_REF].

presented in Section 2. The labels used in column "OEIS" are taken from Sloane's On-Line Encyclopedia of Integer Sequences [START_REF] Sloane | The on-line encyclopedia of integer sequences[END_REF]. The columns "LDE size", resp. "Rec size", refer to the minimal-order homogeneous linear differential, resp. recurrence, equation satisfied by F(t; 1, 1); they contain the order of the equation, and the maximum degree of its polynomial coefficients. The "Pol size" column refers to the algebraicity or transcendence of F(t; 1, 1): cases marked "-" were conjectured transcendental, the other cases were conjectured algebraic and the bidegree of the minimal polynomial was displayed. For example, the generating function F(t; 1, 1) for Kreweras walks (A151265) satisfies a differential equation of order 4 with polynomial coefficients of degree 9 and an algebraic equation P(F(t; 1, 1), t) = 0 for a polynomial P(T, t) of degree 6 in T and 8 in t. The coefficient sequence of F(t; 1, 1) satisfies a recurrence equation of order 6 with polynomial coefficients of degree 4.

For cases 1-22, these conjectural results on D-finiteness, resp. algebraicity, were confirmed by human proofs ‡ obtained almost simultaneously with [START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF] by Bousquet-Mélou and Mishna [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF], using an uniform approach that we will present in §3. We discussed the difficult case 23 (Gessel's model) in §1.15 and §1.16. Concerning the conjectural transcendence results, the first unified proof was given in [START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF] and it is computer-driven; this will be discussed in §1.20. The reference [START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF] also contains the first proof, again computer-driven, that the (differential / recurrence / algebraic) equations conjectured in [START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF] are indeed correct.

As a complement to the results contained in Fig. 7, Bostan and Kauers demon- ‡ Apart from Kreweras' and Gessel's models 20 and 23, the D-finiteness of F S (t; x, y) also follows from: Theorem 6 for the symmetric models 1-16; the Gessel-Zeilberger formula [START_REF] Gessel | Random walk in a Weyl chamber[END_REF] for the "Weyl chamber models" [START_REF] Arquès | Dénombrements de chemins dans R 2 soumis à contraintes[END_REF][START_REF] Askey | How can mathematicians and mathematical historians help each other?[END_REF][START_REF] Ayyer | Towards a human proof of Gessel's conjecture[END_REF][START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF]Th. 2.4] for the "reverse Kreweras model" 21. For the "doubly Kreweras model" 22, [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF]Prop. 15] seems to contain the first proof of D-finiteness, and even of algebraicity. strated that Computer Algebra tools are also able to produce conjectural expressions for the asymptotics of f n = [t n ]F(t; 1, 1). Their results are displayed in Fig. 8 and have been obtained using a combination of algorithmic tools, including Hermite-Padé approximation, constant recognition algorithms built on integer relation detection algorithms like LLL [START_REF] Lenstra | Factoring polynomials with rational coefficients[END_REF] and PSLQ [START_REF] Ferguson | Analysis of PSLQ, an integer relation finding algorithm[END_REF], and convergence acceleration techniques [START_REF] Brezinski | Accélération de la convergence en analyse numérique[END_REF][START_REF] Brezinski | Convergence acceleration during the 20th century[END_REF]. These results have been confirmed a few years later by human proofs by Melczer and Wilson [START_REF] Melczer | Asymptotics of lattice walks via analytic combinatorics in several variables[END_REF], using the theory of analytic combinatorics in several variables [START_REF] Pemantle | Analytic combinatorics in several variables[END_REF]. (Partial results had been previously obtained by Fayolle and Raschel [START_REF] Fayolle | Some exact asymptotics in the counting of walks in the quarter plane[END_REF], Johnson, Mishna and Yeats [START_REF] Johnson | A combinatorial understanding of lattice path asymptotics[END_REF], Duraj [START_REF] Duraj | Random walks in cones: the case of nonzero drift[END_REF], Melczer and Mishna [START_REF] Melczer | Enumerating lattice paths with symmetries through multivariate diagonals[END_REF], Garbit and Raschel [START_REF] Garbit | On the exit time from a cone for random walks with drift[END_REF]).

1.18.

The group of a model. In order to formulate more results on the classification of lattice walks in the quarter plane, we need to introduce an important concept, the group of the walk. To a small-step walk model S one attaches the generating polynomial (also called the inventory) χ S (x, y) := ∑ (i,j)∈S x i y j . This is a bivariate Laurent polynomial in Q[x, x -1 , y, y -1 ], that can be decomposed along powers of x, resp. of y, as follows:

χ S = ∑ (i,j)∈S x i y j = 1 ∑ i=-1 B i (y)x i = 1 ∑ j=-1 A j (x)y j .
The basic, yet fundamental, observation is that χ S (x, y) is left invariant under two rational transformations

ψ(x, y) = x, A -1 (x) A +1 (x) 1 y , φ(x, y) = B -1 (y) B +1 (y) 1 x , y ,
and thus under any element of the group G S := ψ, φ of birational transformations generated by ψ and φ. When it is finite, G S is isomorphic to a dihedral group, since ψ and φ are involutions. This notion of group of a walk originates from a similar notion, introduced in a probabilistic context by Malyshev in the 1970s [START_REF] Malyshev | Wiener-Hopf equations in the quarter-plane, discrete groups and automorphic functions[END_REF]. It was first formally imported in the combinatorial framework by Mishna [START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF][START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF], who realized that the method used in one of Bousquet-Mélou's solutions of the Kreweras model [96, §2.3], the algebraic kernel method, can be used to solve all models with cardinality at most 3. This method is a variation of the classical kernel method: instead of canceling the kernel, it finds a group of actions which fixes the kernel, and which is then used to generate more functional equations that are finally combined together using an algebraic method similar to the reflection principle. Mishna [START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF][START_REF] Mishna | Classifying lattice walks restricted to the quarter plane[END_REF] showed that in the 23 models in Fig. 7, the group is finite, and she determined explicitly its cardinality, which appears to be either 4 (for models 1-16 with an axial symmetry), or 6 (for the models 17, 18, 20, 21, 22, with a diagonal or an anti-diagonal symmetry), or 8 (for the remaining models 19 and 23), see Fig. 9. In a subsequent joint paper, Bousquet-Mélou and Mishna [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF] exploited this idea and managed to solve 22 out of the 23 models in Fig. 7. Their solution will be explained in §3.4. Bousquet-Mélou and Mishna [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF] proved in addition that for all the other 56 models, the group is infinite. Let us sketch their argument, since it is simple, beautiful and very similar to the one used in §1.21. It reduces the question of the infinitude of the group to a (non-)cyclotomy question. Similarly, the argument in §1.21 will reduce the question of non-D-finiteness to the same (non-)cyclotomy question for the same polynomials. (This coincidence, which apparently has not been noticed before, is not fortuitous, see §1.21.) The argument goes as follows. Assume that G S is finite. Then, denoting by θ the composition ψ • φ, the order of θ is finite. Using a Taylor expansion, it follows that for any point (a, b) ∈ C 2 fixed by θ, the order of the Jacobian matrix Jac(θ) at (a, b) is finite, and in particular its two eigenvalues are roots of unity. Now, for all models with infinite group § , there exists a fixed point of θ, and a multiple in Q[t] of the characteristic polynomial of Jac(θ) at that fixed point, that does not contain any cyclotomic factor. This proves that G S is infinite.

At this point, we know that the finiteness of the group for some model implies the D-finiteness of the generating function for that model. One important remaining question is: is the converse true? Another important pending question is: in the Dfinite cases, are there any closed-form expressions for the generating functions? The next two subsections will bring answers and completely clarify the situation. We now focus on models 1-19. The natural question is whether closed-form expressions also exist in these cases. This question has been recently answered in a positive way using Computer Algebra tools in [START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF]: F S is uniformly expressible using iterated integrals of hypergeometric 2 F 1 expressions. More precisely, the following structure result, already conjectured in [84, §3.2], holds true. Note that a similar expression also appears in a related combinatorial context [START_REF] Bostan | Explicit formula for the generating series of diagonal 3D rook paths[END_REF] for rook paths on a three-dimensional chessboard, see Theorem 35 in §3.1.2. Theorem 14 ([76]). Let S be one of the models 1-19 in Fig. 7. Then F S (t; x, y) is expressible as a finite sum of iterated integrals of products of algebraic functions in x, y, t and of expressions of the form 2 F 1 a b c w(t) , where c ∈ N and w(t) ∈ Q(t).

Main results (IV): explicit expressions for models 1-19. Models 20-23 in

Once again, the discovery and the proof of this result are computer-driven; no human proof is available yet. The proof is based, among other tools, on creative telescoping, an efficient algorithmic technique for the symbolic integration of multivariate functions. Details will be discussed in §3.

S occurring 2 F 1 w S occurring 2 F 1 w 1 2 F 1 1 2 , w 16t 2 11 2 F 1 1 2 , 1 2 1 w 16t 2 4t 2 +1 2 2 F 1 1 2 , w 16t 2 12 2 F 1 1 4 , 3 4 1 w 64t 3 (2t+1) (8t 2 -1) 2 3 2 F 1 1 4 , w 64t 2 (12t 2 +1) 2 13 2 F 1 1 4 , 3 4 1 w 64t 2 (t 2 +1) (16t 2 +1) 2 4 2 F 1 1 2 , w 16t(t+1) (4t+1) 2 14 2 F 1 1 4 , 3 4 1 w 64t 2 (t 2 +t+1) (12t 2 +1) 2 5 2 F 1 1 4 , w 64t 4 15 2 F 1 1 4 , 3 4 1 w 64t 4 6 2 F 1 1 4 , w 64t 3 (t+1) (1-4t 2 ) 2 16 2 F 1 1 4 , 3 4 1 w 64t 3 (t+1) (1-4t 2 ) 2 7 2 F 1 1 2 , w 16t 2 4t 2 +1 17 2 F 1 1 3 , 2 3 1 w 27t 3 8 2 F 1 1 4 , w 64t 3 (2t+1) (8t 2 -1) 2 18 2 F 1 1 3 , 2 3 1 w 27t 2 (2t + 1) 9 2 F 1 1 4 , w 64t 2 (t 2 +1) (16t 2 +1) 2 19 2 F 1 1 2 , 1 2 1 w 16t 2 10 2 F 1 1 4 , w 64t 2 (t 2 +t+1) (12t 2 +1) 2
Figure 10. Hypergeometric series occurring in explicit expressions for F(t; x, y). The 2 F 1 are given up to contiguity and derivation, that is, up to integer shifts of the parameters.

The parameters a, b, c of the occurring 2 F 1 's as well as the rational functions w(t) are explicitly given in Table 10. The full expressions of the generating functions F(t; 0, 0), F(t; 0, 1), F(t; 1, 0), F(t; 1, 1), F(t; x, 0), F(t; 0, y) and F(t; x, y) are too large to be displayed here, and are available on-line. It turns out by inspection that the involved hypergeometric functions have a very particular form: they are intimately related to elliptic integrals, namely to the complete elliptic integrals of first and second kinds,

K(k) = π/2 0 (1 -k 2 sin 2 θ) -1/2 dθ = π 2 2 F 1 1 2 , 1 2 1 k 2 , E(k) = π/2 0 (1 -k 2 sin 2 θ) 1/2 dθ = π 2 2 F 1 -1 2 , 1 2 1 k 2 .
For instance, for King walks (case 4), the length generating function is equal to

(9) F(t; 1, 1) = 1 t t 0 1 (1 + 4x) 3 • 2 F 1 3 2 3 2 2 16x(1 + x) (1 + 4x) 2 dx.
See §3.4 for a detailed presentation of this example. Alternatively, an expression of F(t; 1, 1) in terms of elliptic integrals is

F(t; 1, 1) = 1 t t 0 1 π(1 + 4x) 2 x(1 + x) • K 4 x(1 + x) 1 + 4x dx.
The relationship to elliptic integrals appears to hold true in a far more general setting. Indeed, taking Theorem 14 as starting point, van Hoeij has checked that for many (more than 100) integer sequences (a n ) n≥0 in the OEIS whose generating function A(t) = ∑ n≥0 a n t n is both D-finite and convergent in a small neighborhood of t = 0, all second-order irreducible factors of the minimal-order linear differential operator annihilating A(t) are solvable either in terms of algebraic functions, or in terms of complete elliptic integrals. This surprisingly general feature, reminiscent of Dwork's conjecture mentioned in [84, §3.2], begs for a combinatorial explanation.

Main results (V): transcendence for models 1-19.

As said before, models 20-23 in Fig. 7 admit full generating functions that are algebraic. What about the full generating function F S (t; x, y), and its combinatorially meaningful specializations F S (t; 0, 0), F S (t; 1, 0), F S (t; 0, 1), F S (t; 1, 1) for the models 1-23? Computer algebra is able to answer this question. Theorem 15 ([76]). Let S be one of the models 1-19 in Fig. 7. Then for any (α, β) ∈ {(0, 0), (1, 0), (0, 1), (1, 1)}, the power series F S (t; α, β) is transcendental, except in the following four cases:

• S = (model 17) and (α, β) = (1, 1), • S = (model 18) and (α, β) ∈ {(1, 0), (0, 1), (1, 1)}. As a consequence, the power series F S (t; x, y), F S (t; x, 0), and F S (t; 0, y) are transcendental for all the 19 models. Additionally, the generating functions of the four algebraic cases are equal to:

• F (t; 1, 1) = 1 2t 2 1 -t -(1 + t)(1 -3t) , • F (t; 1, 1) = 1 8t 2 1 -2t -(1 + 2t)(1 -6t) , • F (t; 1, 0) = F (t; 0, 1) = 1 32t 3 (1 -6t) 3/2 (1 + 2t) 1/2 -4t 2 + 8t -1 .
Again, the proof of Theorem 15 is computer-driven and crucially relies on the use of several modern Computer Algebra algorithms. This will be discussed in §2.4.5. Algebraicity/transcendence proofs were first considered in some isolated cases: for model [START_REF] Apéry | Sur certaines séries entières arithmétiques[END_REF] They actually proved more, namely that F(t 0 ; x, y) is transcendental for each t 0 ∈ (0, #S -1 ], using the approach in [START_REF] Fayolle | Random walks in the quarter plane[END_REF]Chap. 4]. However, this result does not provide any transcendence information about specializations at x, y ∈ {0, 1}. Note that, for all the 19 models, the excursions generating functions F(t; 0, 0) could alternatively be proved transcendental by an argument based on asymptotics, similar to the one in [START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF]: using results from [START_REF] Denisov | Random walks in cones[END_REF], one can show that the coefficient of t 12n in F(t; 0, 0) grows like κρ n n α for α ∈ {-3, -4, -5}, and this implies transcendence of F(t; 0, 0) by [START_REF] Flajolet | Analytic models and ambiguity of context-free languages[END_REF]Theorem D]. By contrast, note that this asymptotic argument is not sufficient to prove the transcendence of all the other transcendental specializations, as showed for instance by Fig. 8 in the case of F(t; 1, 1) for models 5-10, for which α = -1/2 is not incompatible with algebraicity.

Main results (VI): non-D-finiteness for models with an infinite group.

The last question in view of the complete classification of small step walks in the quarter plane concerns the 56 models with an infinite group. Among them, 5 models are singular; for them, a variant of the kernel method, called the iterated kernel method was used by Mishna and Rechnitzer [START_REF] Mishna | Two non-holonomic lattice walks in the quarter plane[END_REF] (for two models) and by Melczer and Mishna [START_REF] Melczer | Singularity analysis via the iterated kernel method[END_REF] (for all five models), who showed that the length generating function F(t; 1, 1), and thus also the full generating function F(t; x, y), are non-D-finite.

The remaining question concerns the 51 non-singular models with an infinite group: is the full generating function (and its specializations) still non-D-finite?

Computer Algebra is able to help proving the following result.

Theorem 16 ([90]). Let S ⊆ {0, ±1} 2 be any of the 51 nonsingular step sets in N 2 with infinite group G S . Then the generating function F S (t; 0, 0) of S-excursions is not D-finite. Equivalently, the excursion sequence ( f n;0,0 ) n≥0 does not satisfy any nontrivial linear recurrence with polynomial coefficients.

In particular, the full generating function F S (t; x, y) is not D-finite in the 51 cases, since D-finiteness is preserved by specialization [START_REF] Lipshitz | D-finite power series[END_REF]. This corollary had been already obtained by Kurkova and Raschel [START_REF] Kurkova | On the functions counting walks with small steps in the quarter plane[END_REF], but the approach in [START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF] is at the same time simpler, and delivers a more accurate information. This new proof only uses asymptotic information about the coefficients of F S (0, 0, t), and arithmetic information about the constrained behavior of the asymptotics of these coefficients when their generating function is D-finite. More precisely, [START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF] first makes explicit consequences of the general results by Denisov and Wachtel [START_REF] Denisov | Random walks in cones[END_REF] in the case of walks in the quarter plane. This analysis implies that, when n tends to infinity, the excursion sequence f n;0,0 behaves like κ • ρ n • n α , where κ = κ(S) > 0 is a real number, ρ = ρ(S) is an algebraic number, and α = α(S) is a real number such that c =cos( π 1+α ) is an algebraic number. More precisely, [START_REF] André | Solution directe du problème résolu par M. Bertrand[END_REF]

ρ := χ(x 0 , y 0 ), c := ∂ 2 χ ∂x∂y ∂ 2 χ ∂x 2 • ∂ 2 χ ∂y 2 (x 0 , y 0 ), α := -1 -π/ arccos(-c),
where (x 0 , y 0 ) is the unique solution in R 2 >0 of the system ∂χ ∂x

= ∂χ ∂y = 0.
Starting from the step set S, explicit real approximations for ρ, α and c can be determined to arbitrary precision. Moreover, exact minimal polynomials of ρ and c can be determined algorithmically, using tools from elimination theory, namely Gröbner bases [START_REF] Cox | Using algebraic geometry[END_REF]. A classical result in the arithmetic theory of linear differential equations [START_REF] Dwork | An introduction to G-functions[END_REF][START_REF] André | Séries Gevrey de type arithmétique. I. Théorèmes de pureté et de dualité[END_REF][START_REF] Garoufalidis | G-functions and multisum versus holonomic sequences[END_REF] about the possible asymptotic behavior of an integer-valued, exponentially bounded D-finite sequence, states that if such a sequence grows like κ • ρ n • n α , then α is necessarily a rational number. For the 51 cases of nonsingular walks with infinite group, [START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF] proves that the constant α = α(S) is not a rational number. The proof amounts to checking that some explicit polynomials in Q[t] are not cyclotomic. This mirrors the proof of the infinitude of groups for the 51 models, sketched at the end of §1.18. The resemblance is not accidental: with the notations of §1.18, it is possible to prove that (x 0 , y 0 ) is a fixed point for θ and that the characteristic polynomial of the Jacobian Jac(θ) at (x 0 , y 0 ) is equal to T 2 + (2 -4c 2 )T + 1, which admits roots that are roots of unity if and only if α = -1π/ arccos(-c) is a rational number.

E d d T c d d ' E d d T c d d E d d s T
Example 17. Consider the three scarecrows models depicted in Fig. 11. For the first and the third, the approach sketched above shows that the excursions sequence [t n ] F S (t; 0, 0) 1, 0, 0, 2, 4, 8, 28, 108, 372, . . .

is asymptotically equivalent to κ • 5 n • n α , for α = -1 -π/ arccos( 1 4 ) = -3.383396 . . . The irrationality of α prevents F S (t; 0, 0) from being D-finite.
Let us note that a new line of research is currently under development: using a method based on Tutte invariants, Bernardi, Bousquet-Mélou and Raschel [START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF][START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF] showed that for 9 of these 51 models, the generating function is nevertheless Dalgebraic, i.e., it satisfies a system of polynomial (non-linear) differential equations. These models are represented in Fig. 12. In parallel, using differential Galois theory, Dreyfus, Hardouin, Roques and Singer [START_REF] Dreyfus | On the nature of the generating series of walks in the quarter plane[END_REF] proved the hypertranscendence of the remaining 42 models.

1.22. Summary: Classification of 2D non-singular walks. By combining the previous results, we obtain the following classification theorem, which provides a complete characterization of the nonsingular small-step sets with D-finite generating function. Before stating the result, we introduce the notion of orbit sum, that will emerge in §3 in relation with the kernel method. Definition 18. The orbit sum of a quarter-plane model S with finite group G S is the following polynomial in Q[x, x -1 , y, y -1 ]:

OS S := ∑ g∈G S (-1) g g(x)g(y),
where for g ∈ G S we denote by (-1) g the sign of g, which is 1 if g is the product of an even number of generators φ and ψ, and -1 otherwise.

For example in the case of the simple walk OS = x • y -

1 x • y + 1 x • 1 y -x • 1 y .
A simple computation shows that for exactly the four models 20-23, the orbit sum is zero. E.g., for the Kreweras model:

OS = x • y - 1 xy • y + 1 xy • x -y • x + y • 1 xy -x • 1 xy = 0.
We now state the main result of this memoir. Recall that the drift of a model S is defined as the sum of the vectors in S.

Theorem 19. Let S ⊆ {0, ±1} 2 be any of the 74 nonsingular quarter-plane models in Fig. 2. The following assertions are equivalent:

(1) The full generating function F S (t; x, y) is D-finite;

(2) the excursions generating function F S (t; 0, 0) is D-finite;

(3) the excursions sequence

[t 2n ] F S (t; 0, 0) is ∼ K • ρ n • n α , with α ∈ Q; (4
) the group G S is finite;

(5) S has either an axial symmetry, or zero drift and cardinality different from 5.

Moreover, under (1)-( 5), the cardinality of G S is equal to 1)-( 5), F S (t; x, y) is algebraic if and only if S has positive covariance

2 • min ∈ N | α+1 ∈ Z . Still under (
∑ (i,j)∈S ij -∑ (i,j)∈S i • ∑ (i,j)∈S
j > 0 and if and only if OS S = 0. In this case, F S (t; x, y) is expressible using nested radicals.

Otherwise, F S (t; x, y) is expressible using iterated integrals of 2 F 1 expressions.

Proof. Implication (1) ⇒ (2) is easy;

(2) ⇒ (3) is highly non-trivial and follows the combination of a strong probabilistic result [START_REF] Denisov | Random walks in cones[END_REF] and of a strong arithmetic result [START_REF] Dwork | An introduction to G-functions[END_REF][START_REF] André | Séries Gevrey de type arithmétique. I. Théorèmes de pureté et de dualité[END_REF][START_REF] Garoufalidis | G-functions and multisum versus holonomic sequences[END_REF]; (3) ⇒ (4) is the core of the results in [START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF] discussed in §1.21; (4) ⇒ (1) is a consequence of results in [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF][START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF]. The equivalence of ( 2) and ( 5) is read off the tables in Appendix A of [START_REF] Bostan | Non-D-finite excursions in the quarter plane[END_REF]. Condition (5) might seem unnatural; its purpose is to eliminate the three rotations of the "scarecrow" model with step sets depicted in Fig. 11, which have zero drift and non-D-finite generating functions. Finally, the observation on the cardinality can be checked from the data [101, Tables 123].

The characterization of algebraicity in terms of covariance and drift follows by inspection using Theorem 15. The last assertion is Theorem 14.

The classification of walks with small steps in the quarter plane can then be summarized pictorially as follows:

quadrant models S: 79 We conclude this first part of the document with some generalizations and some problems for future investigation.

Walks with unit steps in N 2 . Although small step walks in the quarter plane are quite well understood by now, there remain some open problems. For example, it is still unknown whether the length generating function F(t; 1, 1) is non-D-finite for all 56 models with infinite group. On the other hand, a unified proof is still lacking for the correspondence finite group ↔ D-finite generating function. Walks with unit steps in N 3 . One direction of research concerns the classification of lattice walks in higher dimension. For the moment, an extensive investigation of the case of small step walks in the octant N 3 has been initiated in [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF]. In this case, the notions of the group of a model and of the orbit sum can be mimicked on the 2D case. The first difficulty is the number of cases: there are 2 3 3 -1 ≈ 67 millions models, of which 11 074 225 models are inherently 3-dimensional (instead of 79 in dimension 2). The article [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF] focuses on the 20 804 models that have at most six steps. Among them, 170 cases appear to have a finite group; in the remaining cases, experiments suggest that the group is infinite. Needless to add, Computer Algebra was of crucial help in this study. The full generating function has been proved Dfinite in all the 170 cases, with the exception of 19 intriguing models for which the nature of the generating function still remains unclear. One of them (the fifth in the list below) is the 3D analogue of the Kreweras model.

This leaves an open question: are there 3D non-D-finite models with a finite group? If so, this would constitute a major difference with the 2D case. We have played with the 3D Kreweras model and we conjecture that its generating function is indeed non-D-finite. This is supported by the fact that two different computations suggest that the asymptotics of the sequence k 4n of 3D Kreweras excursions of length 4n (which starts 1, 6, 288, 24444, 2738592, 361998432, . . . ) grows like k 4n ≈ C • 256 n /n 3.3257570041744... , for some C > 0, and the exponent 3.3257570041744 . . . does not appear to be a simple rational number.

Another difference with the case of quarter-plane walks is the disappearance of algebraic models. Certain models do admit algebraic specializations, but then the walks counted by these series do not use all steps of the model, and deleting the unused steps leaves a model of lower dimension. We conjecture that, apart from these degenerate cases, there is no algebraic series among the 3D octant models.

The study [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF] can be summarized as follows.

3D octant models S with ≤ 6 steps: 20 804

|G

non-D-finite?

These results have been recently extended in a computational tour de force by Bacher, Kauers and Yatchak [START_REF] Bacher | Continued classification of 3D lattice models in the positive octant[END_REF] to all 3D octant models: they have found 170 models with |G S | < ∞ and orbit sum 0 (instead of 19 models found by [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF]). Kauers and Wang [START_REF] Kauers | Lattice walks in the octant with infinite associated groups[END_REF] have determined the structure of the group of the models in all these cases, extending results previously obtained by Du, Hou and Wang [START_REF] Du | Infinite orders and non-D-finite property of 3-dimensional lattice walks[END_REF].

Walks with weighted small steps in N 2 . Another line of research concerns the classification of nearest neighbor walks in the quarter plane for models in which Model A Model B Figure 13. Two interesting quadrant models with repeated steps. Both are D-finite, and model B is even algebraic. Note that with only one copy of the repeated step, none of these models would be D-finite ( §1.21).

multiplicities are attached to each direction in the step set. The study has been initiated by Bostan, Bousquet-Mélou, Kauers and Melczer [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF] during their classification of octant models, as it turns out that some 3D models can be reduced by projection to 2D models with multiplicities. Among the octant models, they have identified 14 744 two-dimensional models with at most 6 steps, which yield by projection 527 distinct quadrant models with at most 6 (possibly repeated) steps. Among them, 118 models appeared to have a finite group, of which 95 have a non-zero orbit sum. For 94 of them, the kernel method establishes the D-finiteness of the full generating function, but for one of them (Model A in Fig. 13) Computer Algebra was needed. All the remaining 23 models with finite group and zero orbit sum have been proved algebraic. Among them, 22 can be reduced to a usual quarter-plane model with algebraic generating function, but for the last of them (Model B in Fig. 13) Computer Algebra was needed again. In some sense, models A and B in Fig. 13 are similar to the Gessel model, but much more difficult.

The study in [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF] has been continued by Kauers and Yatchak [START_REF] Kauers | Walks in the quarter plane with multiple steps[END_REF], whose work also heavily relies on Computer Algebra. They carried out a systematic search over all the 4 8 = 65 536 models where each of the eight directions may have any of the four multiplicities 0, 1, 2, 3. Of these, 30 307 were found nontrivial and essentially different. Of these nontrivial models, 1457 turned out to be D-finite (of which 79 models are even algebraic). Of these, three models have a group of order 10, a cardinality that was not possible in the classical (unweighted) setting. Less surprisingly, the correspondence between finite group and D-finite generating function observed in [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF] continues to hold in this weighted 2D context. One open question raised by this study is: does there exist for every n ≥ 2 a quarter-plane model with multiplicities whose group has order 2n? This is true in a probabilistic context, but for a different notion of group [START_REF] Fayolle | Random walks in the quarter plane[END_REF]. In a very recent work, Courtiel, Melczer, Mishna and Raschel [START_REF] Courtiel | Weighted lattice walks and universality classes[END_REF] push even further the investigation of weighted models. Other extensions. There are many other questions on the combinatorics of lattice paths in the cones, and certainly Computer Algebra will have a word to say, at least for some of them. Counting walks in non-convex cones is currently under investigation: after the case of the slit plane [START_REF] Bousquet-Mélou | Counting paths on the slit plane[END_REF][START_REF] Bousquet-Mélou | Walks on the slit plane: other approaches[END_REF][START_REF] Bousquet-Mélou | Walks on the slit plane[END_REF][START_REF] Xin | Proof of a conjecture on the slit plane problem[END_REF][START_REF] Rubey | Transcendence of generating functions of walks on the slit plane[END_REF], it is now the turn of the cone C := {(i, j) : i ≥ 0 or j ≥ 0} [START_REF] Bousquet-Mélou | Square lattice walks avoiding a quadrant[END_REF]. Also, walks with larger steps in the quadrant are currently under investigation [START_REF] Fayolle | About a possible analytic approach for walks in the quarter plane with arbitrary big jumps[END_REF][START_REF] Bostan | Walks with large steps in an orthant[END_REF]. There are several challenges, among them to find (and to use!) a notion close to the group of a model, which was specific to small-step models. For instance, for S = {(0, 1), (1, -1), (-2, -1)}, Bostan, Bousquet-Mélou and Melczer show that such a notion exists, and allows to prove that F S (t; x, y) is D-finite, via the positive part representation:

xyF S (t; x, y) = [x >0 y >0 ] (x -2x -2 )(y -(x -x -2 )y -1 ) 1 -t(xy -1 + y + x -2 y -1
) .

Guess-and-Prove.

What is "scientific method"? Philosophers and non-philosophers have discussed this question and have not yet finished discussing it. Yet as a first introduction it can be described in three syllables: Guess and test. Mathematicians too follow this advice in their research although they sometimes refuse to confess it. They have, however, something which the other scientists cannot really have. For mathematicians the advice is First guess, then prove.

G. Pólya [START_REF] Pólya | Guessing and proving[END_REF].

In this second part of the document, we enter into more technical details related to the experimental mathematics methodology that was employed to discover and to prove an important part of the results presented in §1, notably related to the celebrated Gessel walks ( §1.15, §1.16) and more generally to the classification of lattice path models with D-finite generating functions ( §1.17, §1.23). The process of experimental mathematics is to discover mathematical phenomenona by observing them via computations, before formally proving them. The rigorous proving step may be human, in the spirit of classical mathematics, or itself computerized, in the spirit of the current memoir. One of the experimental mathematics paradigms that was intensively used in recent years in the lattice path combinatorics context is the so-called guess-and-prove approach. It was introduced in this combinatorial context in work by Bostan and Kauers [START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF][START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF], but its roots can be found in Pólya's remarkable books [START_REF] Pólya | How to solve it[END_REF][START_REF] Pólya | Mathematics and plausible reasoning[END_REF], who popularized it as a fruitful mathematical proof strategy. The power of the method is highly enhanced when used on a computer, in conjunction with fast Computer Algebra algorithms.

This enhancement could be called the automated (or, algorithmic) guess-and-prove approach, and it is the topic of the current section. The first half (the guessing part) of the approach is based on a "functional interpolation" phase, which consists in recovering equations starting from (truncations of) solutions. It is called differential approximation [START_REF] Guttmann | On a new method of series analysis in lattice statistics[END_REF][START_REF] Khan | High-order differential approximants[END_REF], or algebraic approximation [START_REF] Brak | Algebraic approximants: a new method of series analysis[END_REF], depending on the type of equations to be reconstructed. For instance, differential approximation is an operation to get an ODE likely to be satisfied by a given approximate series expansion of an unknown function. This has been used at least since the 1970s by physicists [START_REF] Guttmann | On a new method of series analysis in lattice statistics[END_REF][START_REF] Guttmann | On the recurrence relation method of series analysis[END_REF], and is a key stone in recent spectacular applications in experimental mathematics, such as [START_REF] Koutschan | Proof of George Andrews's and David Robbins's q-TSPP conjecture[END_REF]. Modern versions [START_REF] Salvy | Gfun: a Maple package for the manipulation of generating and holonomic functions in one variable[END_REF][START_REF] Kauers | Guessing handbook[END_REF][START_REF] Hebisch | Extended rate, more GFUN[END_REF] are based on subtle algorithms for Hermite-Padé approximants [START_REF] Beckermann | A uniform approach for the fast computation of matrix-type Padé approximants[END_REF]. The second half (the proving part) of the approach is based on fast manipulations (e.g., resultants and factorization) with exact algebraic objects (e.g., polynomials and differential operators).

Methodology for proving algebraicity and D-finiteness.

We illustrate the general principles of the guess-and-prove method when applied to proving that, for some lattice path model S with small steps in the quarter plane, the full generating function F S (t; x, y) is D-finite or algebraic. Recall from §1.14 that the problem amounts to solving the kernel equation ( 7): K(t; x, y)F(t; x, y) = xy + K(t; x, 0)F(t; x, 0) + K(t; 0, y)F(t; 0, y) -K(t; 0, 0)F(t; 0, 0), where K S (t; x, y) = xy(1 -t∑ (i,j)∈S x i y j ) is the kernel polynomial.

The method can be decomposed into three main steps: (S1) Data generation: one first computes a high order expansion of the power series F S (t; x, y);

(S2) Conjecture: from the local information computed at Step (S1), one tries to guess a global information, namely a candidate for a polynomial, resp. for a system of linear differential equations, satisfied by F S (t; x, y); this is done by using algebraic, resp. differential, approximation; (S3) Proof: one rigorously certifies the output of Step (S2), by using (exact) computations on multivariate polynomials, and on linear differential equations with polynomial coefficients. In practice, Steps (S1), (S2), (S3) are performed using efficient algorithms from Computer Algebra.

As it turns out, an important improvement from the complexity of computations viewpoint is to perform the guessing step (S2) on the sections F S (t; x, 0) and F S (t; 0, y) only. This is sufficient due to the kernel equation, since both algebraicity and D-finiteness are preserved by sums, products and specializations. This simplification is crucial, as equations for the sections are usually much smaller than equations for the full generating function.

In §2.2, §2.3 and §2.4 we take a closer look at Steps (S1), (S2) and (S3).

2.2.

Step (S1): high order series expansions. The first step of the method is based on a very basic observation: the full counting sequence ( f n;i,j ) satisfies a recurrence with constant coefficients (11)

f n+1;i,j = ∑ (k, )∈S
f n;i-k,j-for n, i, j ≥ 0 with the initial conditions f 0;i,j = δ 0,i,j and f n;-1,j = f n;i,-1 = 0. The recurrence simply translates the step-by-step construction of quarter plane walks with model S: a S-walk of length n finishing at (i, j) is obtained from a walk of length n -1, followed by a step in S; the initial conditions translate the quarter-plane constraint.

Notice that as in the case of the much simpler kernel equation ( 2), multiplying the recurrence (11) by t n x i y j , summing over n, i, j, and using the initial conditions yields the kernel equation [START_REF] Albert | Permutations sortable by two stacks in parallel and quarter plane walks[END_REF].

Example 20. For the Kreweras walks, where k n;i,j denotes f n;i,j for S = ,

k n+1;i,j = k n;i+1,j + k n;i,j+1 + k n;i-1,j-1 .
The recurrence [START_REF] André | G-functions and geometry[END_REF] can be used to determine the value of f n;i,j for specific integers n, i, j ∈ N. The inequality f n;i,j ≤ #S n implies that f n;i,j is a non-negative integer whose bit size is at most O(n). Therefore, if N ∈ N, the truncated power series F S (t; x, y) mod t N can be computed by a straightforward algorithm that uses O(N 3 ) arithmetic operations and Õ(N 4 ) bit operations. (We assume that two integers of bit-size N are multiplied in Õ(N) bit operations using FFT [START_REF] Schönhage | Schnelle Multiplikation großer Zahlen[END_REF]; here, the soft-O notation Õ( ) hides logarithmic factors.) The memory storage requirement is proportional to N 3 . For the experiments made in [START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF], N = 1000 was chosen. With this choice, the computation of the f n;i,j becomes very time and memory consuming.

Example 21. For the Kreweras model, one obtains

K(t; x, y) = 1 + xyt + (x 2 y 2 + y + x)t 2 + (x 3 y 3 + 2xy 2 + 2x 2 y + 2)t 3 + (x 4 y 4 + 3x 2 y 3 + 3x 3 y 2 + 2y 2 + 6xy + 2x 2 )t 4
+ (x 5 y 5 + 4x 3 y 4 + 4x 4 y 3 + 5xy 3 + 12x 2 y 2 + 5x 3 y + 8y + 8x)t 5 + • • • , from which the first terms of the length generating function K(t; 1, 1) are computed

K(t; 1, 1) = 1 + t + 3t 2 + 7t 3 + 17t 4 + 47t 5 + 125t 6 + 333t 7 + 939t 8 + 2597t 9 + 7183t 10 + 20505t 11 + 57859t 12 + 163201t 13 + 469795t 14 + • • • .
To summarize, step (S1) is very simple mathematically, but the naive algorithm used for it is not satisfactory. Its weakness is that in order to compute an univariate series such as F S (t; 1, 1), or a bivariate series like F S (t; x, 0), it needs to compute the trivariate series F S (t; x, y). An important problem is to accelerate this algorithm. Our suggestion is to devise a divide-and-conquer method based on eq. ( 18) on p. 38, in the spirit of the algorithms in [START_REF] Brent | Fast algorithms for manipulating formal power series[END_REF][START_REF] Brent | On the complexity of composition and generalized composition of power series[END_REF][START_REF] Bostan | Fast computation of power series solutions of systems of differential equations[END_REF][START_REF] Bostan | Power series solutions of singular (q)-differential equations[END_REF]. This would allow to compute the sections F S (t; x, 0) mod t N and F S (t; 0, y) mod t N in quasi-optimal time (i.e., almost linear in their size, up to logarithmic factors), from which F S (t; 1, 1) mod t N could be easily reconstructed using the kernel equation ( 7) evaluated at x = y = 1.

2.3.

Step (S2): guessing equations. The purpose of the second step of the method is to guess (differential, or algebraic) equations for F S (t; x, y).

A first idea.

A first, but crucial, simplification comes from the simple remark that the kernel equation [START_REF] Albert | Permutations sortable by two stacks in parallel and quarter plane walks[END_REF] expresses the full generating function F S (t; x, y) as a linear combination with rational function coefficients in Q(x, y, t) of its sections F S (t; x, 0), F S (t; 0, y) and F S (t; 0, 0). Therefore, by closure properties of algebraic and D-finite functions [START_REF] Lipshitz | D-finite power series[END_REF], F S (t; x, y) is D-finite (resp., algebraic) if and only if its sections F S (t; 0, y) and F S (t; 0, 0) are both D-finite (resp., algebraic).

Example 22. In terms of generating functions, the recurrence in Ex. 20 reads [START_REF] André | Séries Gevrey de type arithmétique. I. Théorèmes de pureté et de dualité[END_REF] xy -(x + y + x 2 y 2 )t K(t; x, y) = xyxt K(t; x, 0)yt K(t; 0, y).

In order to prove the D-finiteness, resp. the algebraicity, of K(t; x, y), it is enough to prove the D-finiteness, resp. the algebraicity, of its sections K(t; x, 0) and K(t; 0, y).

In some cases, this simplification is crucial; for instance, in the case of the Gessel model, the minimal polynomial of F(t; x, y) has a size of ≈ 30Gb, while sizes of the minimal polynomials of the sections F(t; x, 0) and F(t; 0, y) are merely ≈ 1Mb.

Guessing equations for the sections F S (t;

x, 0) and F S (t; 0, y). At the end of Step (S1), we are reduced to performing the following guessing tasks. Task 1 (differential guessing): Given the first N terms of S = F S (t;

x, 0) ∈ Q[x][[t]],
search for a linear differential equation satisfied by S at precision N:

(13) L x,0 (S) = c r (x, t) • ∂ r S ∂t r + • • • + c 1 (x, t) • ∂S ∂t + c 0 (x, t) • S = 0 mod t N .
Task 2 (algebraic guessing): Given the first N terms of S = F S (t;

x, 0) ∈ Q[x][[t]],
search for a polynomial equation satisfied by S at precision N:

(14) P x,0 (S) = c r (x, t) • S r + • • • + c 1 (x, t) • S + c 0 (x, t) • 1 = 0 mod t N .
Here and below, we use the compact notation P x,0 for a trivariate polynomial in Q[T, t, x], and L x,0 for an operator in the Weyl algebra Q(x, t) ∂ t of linear differential operators in ∂ t = d dt with rational function coefficients in Q(t, x). We use the similar notation L 0,y (S ) and P 0,y (S ) for equations potentially satisfied by the other section S = F S (t; 0, y)

∈ Q[y][[t]].
The idea behind differential guessing is that if the given power series S happens to be D-finite, then for a sufficiently large N, a differential equation of type [START_REF] André | Sur la conjecture des p-courbures de Grothendieck-Katz et un problème de Dwork[END_REF] (thus satisfied a priori only at precision N) will provide a differential equation which is really satisfied by S in Q[x][[t]] (i.e., at precision infinity). In other words, the (conjectural) D-finiteness of a power series can be eventually recognized using a finite amount of information. The same holds for the algebraic guessing.

Example 23 (continued). Using N = 80 terms of K(t; x, 0) = F (t; x, 0), one can guess a linear differential operator of order 4, and degrees [START_REF] Apagodu | Multi-variable Zeilberger and Almkvist-Zeilberger algorithms and the sharpening of Wilf-Zeilberger theory[END_REF][START_REF] André | G-functions and geometry[END_REF] in (t, x):

L x,0 = t 3 • (3t -1) • (9t 2 + 3t + 1) • (3t 2 + 24t 2 x 3 -3xt -2x 2 )• • (16t 2 x 5 + 4x 4 -72t 4 x 3 -18x 3 t + 5t 2 x 2 + 18xt 3 -9t 4 )• • (4t 2 x 3 -t 2 + 2xt -x 2 ) • ∂ 4 t + • • • such that L x,0 (K(t; x, 0)) = 0 mod t 80 .
Similarly, one can guess a polynomial of degree (6, 10, 6) in (T, t, x)

P x,0 = x 6 t 10 T 6 -3x 4 t 8 (x -2t)T 5 + + x 2 t 6 12t 2 + 3t 2 x 3 -12xt + 7 2 x 2 T 4 + • • • such that P x,0 (K(t; x, 0)) = 0 mod t 80 .
Therefore, it is very likely that K(t; x, 0) verifies the linear differential equation L x,0 (K(t; x, 0)) = 0 and the algebraic equation P x,0 (K(t; x, 0)) = 0, but at this stage we only have experimental evidence, which is by no means a rigorous proof.

In Tasks 1 and 2, the unknowns c i are (not simultaneously zero) polynomials in Q[x, t]. If their degrees in t are bounded by some prescribed integer d ≥ 0 such that (d + 1)(r + 1) > N, then a simple linear algebra argument shows that a differential equation of type [START_REF] André | Sur la conjecture des p-courbures de Grothendieck-Katz et un problème de Dwork[END_REF], resp. an algebraic equation of type [START_REF] Apagodu | Multi-variable Zeilberger and Almkvist-Zeilberger algorithms and the sharpening of Wilf-Zeilberger theory[END_REF], should exist. On the other hand, if d, r and N are such that (d + 1)(r + 1) N, then equation ( 13) and ( 14) translate into highly over-determined linear systems, which have no reason to possess non-trivial solutions, unless S really is D-finite, resp. algebraic.

All previous remarks also apply to any specialization of S to same value x ∈ Q. The pending question is: how to solve efficiently Tasks 1 and 2, given d, r, N? Obviously, both amount to solving linear algebra problems in size N over Q(x). More precisely, a candidate differential, resp. polynomial, equation of type [START_REF] André | Sur la conjecture des p-courbures de Grothendieck-Katz et un problème de Dwork[END_REF], resp. ( 14), for S can be computed by Gaussian elimination. But the corresponding systems are not randomly dense linear systems. They possess a very special structure, that can be exploited algorithmically in several ways. First, instead of solving linear systems of size N over Q(x), it is better to use an evaluation-interpolation scheme: evaluate the system at several points x, solve the corresponding systems over Q, and recombine the results by interpolation. The evaluation and interpolation steps can be performed very efficiently [START_REF] Zur Gathen | Modern computer algebra[END_REF]Chap. 10], especially at points in geometric progression [91, §5]. Second, instead of solving linear systems over Q, it is better to solve several systems over finite fields F p using a modular approach: the linear algebra step is performed modulo several primes p, and the results are recombined over Q via rational reconstruction based on an effective version of the Chinese remainder theorem. Again, this can be performed very efficiently [START_REF] Zur Gathen | Modern computer algebra[END_REF]Chap. 10]. Third, instead of using Gaussian elimination for solving the linear systems over F p that arise from ( 13) and ( 14) by specialization and reduction, it is better to exploit their Toeplitzlike structure: their matrices are obtained by concatenation of Sylvester-like blocks, that possess the Toeplitz property of diagonal invariance, see §2.5 for details. Said differently, equations ( 13) and ( 14) are particular instances of Hermite-Padé approximation problems, and can be solved very efficiently. More precisely, while Gaussian elimination in size N over F p has cubic arithmetic complexity in N, fast algorithms for Hermite-Padé approximation have quasi-linear complexity in N, see §2.5.3. Such sophisticated algorithms rely on fast (FFT-like) arithmetic for the polynomial ring F p [t] [START_REF] Zur Gathen | Modern computer algebra[END_REF][START_REF] Bini | Polynomial and matrix computations[END_REF][START_REF] Bürgisser | Algebraic complexity theory[END_REF] In practical implementations, for a given precision N, one searches for equations of increasing order r = 1, 2, . . ., and a corresponding degree d ≈ N/r. If no differential equation like ( 13) is found, this definitely rules out the possibility that a differential equation of order r and degree d exists. However, this does not imply that the series at hand is not D-finite. It may still be that S satisfies a differential equation of order higher than r, or an equation with polynomial coefficients of degree exceeding d. In that case, one doubles the series precision N, and starts over.

Sometimes (see §2.3.3 and §2.4.5) one needs to obtain the minimal-order differential equation L min (S) = 0 satisfied by the given generating power series S. The choice (d, r) of the target degree and order does not necessarily lead to the minimal operator L min . Worse, it may even happen that the number of initial terms N is not large enough to allow the recovery of L min , while these N terms suffice to guess non-minimal order operators. The explanation of why such a situation occurs systematically was first given in [START_REF] Bostan | Differential equations for algebraic functions[END_REF] for the case of differential equations satisfied by algebraic functions: minimal-order differential equations are often cluttered with apparent singularities, which considerably increase the degree of their coefficients. Therefore, they require too many terms N of the series S, and this prevents, or slows down, the reconstruction of equations. Differential guessing can benefit from the calculation of non-minimal equations, by minimizing not the order but the total size of the output. These considerations are intimately related to the operation of desingularization [START_REF] Chen | Order-degree curves for hypergeometric creative telescoping[END_REF][START_REF] Chen | Trading order for degree in creative telescoping[END_REF][START_REF] Chen | Desingularization explains order-degree curves for Ore operators[END_REF][START_REF] Chen | Desingularization of Ore operators[END_REF]. All in all, a good heuristic to get L min is to compute several non-minimal operators and to take their greatest common right divisor (gcrd); generically, the result is exactly L min .

Example 24. For Gessel walks, N = 1000 terms of G(t; x, y) = F (t; x, y) are sufficient to guess

• a differential operator L x,0 ∈ Q(x, t) ∂ t , of order 11 in ∂ t , bidegree [START_REF] Bousquet-Mélou | Walks in the quarter plane: Kreweras' algebraic model[END_REF][START_REF] Bostan | Algebraic diagonals and walks[END_REF] in (t, x), and integer coefficients of at most 61 digits • a differential operator L 0,y ∈ Q(y, t) ∂ t , of order 11 in ∂ t , bidegree [START_REF] Bostan | Fast computation of the Nth term of an algebraic series over a finite prime field[END_REF][START_REF] Barton | Some aspects of the random sequence[END_REF] in (t, y), and integer coefficients of at most 51 digits such that L x,0 (G(t; x, 0)) = L 0,y (G(t; 0, y)) = 0 mod t 1000 .

Here is the way this was done. For a fixed value a, and modulo a fixed prime p, several (non-minimal order) operators in F p [t] ∂ t for G(t; a, 0) can be guessed by Guessing differential operators for G(t; a, 0), for prime p and a ∈ F p : minimal-order operator (blue point above the hyperbola) obtained as gcrd of several non-minimal operators (blue points below the hyperbola). Points below the hyperbola correspond to operators obtainable by Hermite-Padé approximation with 1000 terms.

Hermite-Padé approximation using 1000 terms of G(t; a, 0). Some of them are represented by the blue points below the hyperbola in Figure 14, e.g., one of them has order 14 and degree 43 in t. However, interpolating from one of those an operator in Q[t, x] ∂ t for G(t; x, 0) appears to be computationally extremely expensive. The reconstruction (w.r.t. x) becomes feasible (in reasonable degree 78) when applied to the minimal-order operators (represented by the blue point above the hyperbola), themselves obtained as gcrds in F p (t) ∂ t of several non-minimal operators. Note that without gcrds, the minimal-order operator could not have been found by Hermite-Padé approximation with only 1000 terms. Also note that guessing L x,0 naively by undetermined coefficients would have required solving a dense linear system of size 91956 with ≈ 1000 digits entries! As a historical note, the discovery in 2008 of L x,0 and L 0,y first led Bostan and Kauers [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF] suspect that G(t; x, y) is D-finite.

Efficient implementation of differential and algebraic guessing procedures have been implemented in most computer algebra systems, see e.g., the Maple package gfun written by Salvy and Zimmermann [START_REF] Salvy | Gfun: a Maple package for the manipulation of generating and holonomic functions in one variable[END_REF], the Mathematica package Guess.m by Kauers [START_REF] Kauers | Guessing handbook[END_REF], or the FriCAS package Guess written by Hebisch and Rubey [START_REF] Hebisch | Extended rate, more GFUN[END_REF].

Empirical certification of guesses.

Confidence in guessed equations can be complemented by using various filters. Once discovered a differential equation [START_REF] André | Sur la conjecture des p-courbures de Grothendieck-Katz et un problème de Dwork[END_REF] or an algebraic equation ( 14) that the power series S seems to satisfy, it is useful to inspect several properties of these equations, in order to provide more convincing evidence that they are correct. These properties have various flavors: algebraic, analytic and even arithmetic. If the candidate guessed equations pass these filters, this offers striking experimental evidence that they are not artefacts.

Algebraic sieve: High order series matching. The equations ( 13) and ( 14) are obtained starting from N coefficients of the power series S. They are therefore satisfied a priori only modulo t N . One can compute more terms of S, say 2N, and check whether the same equations still hold modulo t 2N . If this is the case, chances increase that the guessed equations also hold at infinite precision.

Analytic sieve: Singularity analysis. For any a ∈ N, the univariate power series F S (t; a, 0) has integer coefficients and positive radius of convergence. Thus, if in addition it is D-finite, then it is a G-function [START_REF] Dwork | An introduction to G-functions[END_REF]. General results by Katz and Honda [START_REF] Katz | Nilpotent connections and the monodromy theorem: Applications of a result of Turrittin[END_REF][START_REF] Honda | Algebraic differential equations[END_REF], and Chudnovsky [START_REF] Chudnovsky | Applications of Padé approximations to Diophantine inequalities in values of G-functions[END_REF] then imply that the minimal order differential operator for F S (t; a, 0) needs to be Fuchsian (it admits only regular singularities, including at infinity) and its exponents at each singularity must be rational numbers. See [START_REF] André | G-functions and geometry[END_REF][START_REF] Chambert-Loir | Théorèmes d'algébricité en géométrie diophantienne[END_REF][START_REF] Dwork | An introduction to G-functions[END_REF] for more details on this topic. Arithmetic sieve: G-functions and global nilpotence. Last, but not least, one may 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 [START_REF] Chudnovsky | Applications of Padé approximations to Diophantine inequalities in values of G-functions[END_REF], the minimal order differential operator L = L S min ∈ Q[t] ∂ t killing a G-function S 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 µ ≥ 1 such that the remainder of the Euclidean (right) division of ∂ pµ t by L is congruent to zero modulo p [START_REF] Honda | Algebraic differential equations[END_REF][START_REF] Dwork | Differential operators with nilpotent p-curvature[END_REF]. From a computational view-point, a fine feature is that the nilpotence modulo p is checkable. If r denotes the order of L, let A p (L) be the p-curvature matrix of L, defined as the r × r matrix with entries in Q(t) whose (i, j) entry is the coefficient of ∂ j-1 t in the remainder of the Euclidean (right) division of ∂ p+i-1 t by L. Then, L is nilpotent modulo p if and only if the matrix A p (L) is nilpotent modulo p [START_REF] Dwork | Differential operators with nilpotent p-curvature[END_REF][START_REF] Schmitt | Operators with nilpotent p-curvature[END_REF]. Faster tests exist [START_REF] Bostan | Fast algorithms for differential equations in positive characteristic[END_REF][START_REF] Bostan | A fast algorithm for computing the characteristic polynomial of the p-curvature[END_REF][START_REF] Bostan | A fast algorithm for computing the p-curvature[END_REF][START_REF] Bostan | Computation of the similarity class of the p-curvature[END_REF].

This yields a fast algorithmic filter: as soon as we guess a candidate differential equation satisfied by a generating function which is suspected to be a G-function (e.g., by F(t; 1, 1)), we check whether its p-curvature matrix A p (L) is nilpotent, say modulo the first 50 primes for which the reduced operator L mod p is well-defined. If A p (L) is indeed nilpotent modulo p for all those primes p, then the guessed equation is, with very high probability, the correct one. This arithmetic sieving can be pushed even further. A famous conjecture, attributed to Grothendieck [START_REF] Katz | Algebraic solutions of differential equations (p-curvature and the Hodge filtration)[END_REF][START_REF] Katz | A conjecture in the arithmetic theory of differential equations[END_REF][START_REF] André | Sur la conjecture des p-courbures de Grothendieck-Katz et un problème de Dwork[END_REF], asserts that the differential equation L(S) = 0 possesses a basis of algebraic solutions (over Q(t)) if and only if A p (L) is zero modulo p for almost all primes p. Even if the conjecture is, for the moment, fully proved only in special cases [START_REF] Chambert-Loir | Théorèmes d'algébricité en géométrie diophantienne[END_REF] (notably for Picard-Fuchs equations [START_REF] Katz | Algebraic solutions of differential equations (p-curvature and the Hodge filtration)[END_REF]) one can use it as an oracle to detect whether a guessed differential equation has a basis of algebraic solutions.

Example 25 (continued). For Gessel walks, the guessed (order-11) operators L x,0 and L 0,y for G(t; x, y) = F (t; x, y) pass all the preceding filters, including the one based on p-curvatures. More precisely, for randomly chosen prime number p, and a, b ∈ F p , both L a,0 and L 0,b right-divide the pure power ∂ 11•p t in F p (x) ∂ t . These operators actually have a stronger property: they even right-divide ∂ p t ; in other terms, they have zero p-curvature for all the tested primes p. This was the key observation in the discovery [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF] that the trivariate generating function for Gessel walks is algebraic.

The reader may wonder why the authors of [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF] did not try algebraic guessing first. The first reason is that they had no reason to suspect that G(t; x, y) is algebraic, since even the specialization G(t; 0, 0) was generally thought to be transcendental. The second reason is that more terms of G(t; x, y) are needed to recognize alebraicity (1200, instead of 1000, see below), and the power series expansion to such high orders is computationally very expensive both in time and memory.

Example 26 (continued). Still for Gessel walks, now using N = 1200 terms of G(t; x, y) = F (t; x, y) is sufficient to guess annihilating polynomials for sections:

• P x,0 ∈ Z[T, t, x] of degree [START_REF] Ball | Irrationalité d'une infinité de valeurs de la fonction zêta aux entiers impairs[END_REF][START_REF] Bini | Polynomial and matrix computations[END_REF][START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF], integer coefficients of at most 21 digits, • P 0,y ∈ Z[T, t, y] of degree [START_REF] Ball | Irrationalité d'une infinité de valeurs de la fonction zêta aux entiers impairs[END_REF][START_REF] Blanchard | Random walks and diffusions on graphs and databases[END_REF][START_REF] Beukers | Monodromy for the hypergeometric function n F n-1[END_REF], integer coefficients of at most 23 digits, such that P x,0 (G(t; x, 0)) = P 0,y (G(t; 0, y)) = 0 mod t 1200 .

Step (S3): rigorous proof.

Guessing is good, proving is better.

G. Pólya [START_REF] Pólya | How to solve it[END_REF].

Basic idea.

The third and last step of the guess-and-prove method (for a quarter-plane model S for which the first two steps are assumed to have succeeded) consists in rigorously proving that the candidate (guessed) equations are indeed correct. Roughly, the basic idea is the following. Assume that one guessed equations for F S (t; x, y) which admit a solution S(t; x, y) in some power series ring R, typically

Q[[x, y, t]] or Q[x, x -1 , y, y -1 ][[t]
], in which the kernel equation ( 7) has a unique solution, namely F S (t; x, y). Then, using effective closure properties for algebraic and D-finite functions [START_REF] Lipshitz | D-finite power series[END_REF] enables to compute the same (algebraic, or differential) equations for both sides of the kernel equation ( 7) with F S (t; x, y) replaced by S(t; x, y), and to eventually prove that the identity [START_REF] Apéry | Sur certaines séries entières arithmétiques[END_REF] K(t; x, y)S(t; x, y) = xy + K(t; x, 0)S(t; x, 0) + K(t; 0, y)S(t; 0, y) -K(t; 0, 0)S(t; 0, 0) holds in R, where K S (t; x, y) = xy(1 -t∑ (i,j)∈S x i y j ). By uniqueness, it follows that F S (t; x, y) and S(t; x, y) coincide, and thus F S (t; x, y) is indeed algebraic (or D-finite), since S(t; x, y) is so, by design.

In practice, contrary to this ideal scenario, equations for the full generating function are too big to be computed, at least in many interesting cases. As explained in §2.3, one only has access to guessed equations for the sections F S (t; x, 0) and F S (t; 0, y). In this case, a variant of the method is used, and it is based on the reduced kernel equation, see §2.4.3 below. But before going into this, let us illustrate the guess-and-prove philosophy on a simpler example.

Warm-up: algebraicity of Gessel excursions.

Let us prove that the generating function G(t; 0, 0) of Gessel excursions is algebraic, by taking Theorem 11 as the starting point. In other words, let us prove the algebraicity of the power series

g(t) := G( √ t; 0, 0) = ∞ ∑ n=0 (5/6) n (1/2) n (5/3) n (2) n (16t) n .
Of course, one could appeal to a proof that relies on equation ( 8) and on Schwarz's classification [START_REF] Schwarz | Über diejenigen Fälle, in welchen die Gaußische hypergeometrische Reihe einer algebraische Funktion ihres vierten Elementes darstellt[END_REF] of algebraic 2 F 1 s, or other methods discussed in §1.11, like the Beukers-Heckman criterion (Theorem 8). Compared to these proofs, the constructive proof given below has the advantage that it can be applied similarly in situations where no classification results are available.

The guess-and-proof method works as follows: first guess a polynomial P(t, T) in Q[t, T], then prove that P admits the power series g(t) = ∑ ∞ n=0 g n t n as a root, where g n = (5/6) n (1/2) n

(5/3) n (2) n 16 n . In more details, the proof decomposes into three main steps:

1. (Guessing) A suitable polynomial P (see below) can be guessed automatically from the first 100 terms of g(t) using the approach explained in §2. 3. (Proof) r(t) = ∑ ∞ n=0 r n t n being algebraic, it is D-finite ( §1.11), and thus its coefficients satisfy a recurrence with polynomial coefficients, which is [START_REF] Apéry | Irrationalité de ζ(2) et de ζ(3)[END_REF] (n + 2)(3n + 5)r n+1 -4(6n + 5)(2n + 1)r n = 0, r 0 = 1.

Thus r n = (5/6) n (1/2) n

(5/3) n (2) n 16 n = g n , and g(t) = r(t) is algebraic. The concrete computations can be performed for instance in Maple using the package gfun, which provides the commands algeqtodiffeq for the algebraic guessing task in Step 1, algeqtodiffeq for Abel's theorem in Step 3 and diffeqtorec for the conversion differential equation → recurrence in Step 3. The result of the two lines ¶ > P:=gfun:-listtoalgeq([seq(pochhammer(5/6,n)*pochhammer(1/2,n)/ pochhammer(5/3,n)/pochhammer(2,n)*16^n, n=0..100)], g(t)): > gfun:-diffeqtorec(gfun:-algeqtodiffeq(P [START_REF] Abel | OEuvres complètes. Tome II[END_REF], g(t)), g(t), r(n)); is the recurrence [START_REF] Apéry | Irrationalité de ζ(2) et de ζ(3)[END_REF].

Algebraicity proofs for Kreweras and Gessel walks.

We now sketch the last part of the guess-and-prove method for proving the algebraicity of the generating functions for Kreweras and for Gessel walks. We focus on the Kreweras case, since the computations are easier and most of the ideas are already present.

The proof follows the same principles as the one just explained in §2.4.2. The idea is to guess, then to certify annihilating polynomials. The main difference with the situation in §2.4.2 is that an explicit closed form expression is no longer available beforehand for the power series whose algebraicity needs to be proved. Instead, we only have implicit equations that define that series. The method has three steps, and consists in applying the basic idea explained in §2.4.1, with the major difference that we cannot afford guessing of equations for the full generating function. The first step produces a so-called reduced kernel equation for the sections F(t; x, 0) and F(t; 0, y). In the Kreweras case, the step set being symmetric with respect to the main diagonal, the generating function K(t; x, y) enjoys the property K(t; x, y) = K(t; y, x), which simplifies the kernel equation ( 7) to [START_REF] Arquès | Dénombrements de chemins dans R 2 soumis à contraintes[END_REF] (xy -(x + y + x 2 y 2 )t)K(t; x, y) = xy -xtK(t; x, 0) -ytK(t; y, 0).

The proof goes as follows (the corresponding computations, performed in Maple, are displayed in Fig. 15):

# HIGH ORDER EXPANSION (S1) > st,bu:=time(),kernelopts(bytesused): > f:=proc(n,i,j) option remember; if i<0 or j<0 or n<0 then 0 elif n=0 then if i=0 and j=0 then 1 else 0 fi else f(n-1,i-1,j-1)+f(n-1,i,j+1)+f(n-1,i+1,j) fi end: > S:=series(add(add(f(k,i,0)*x^i,i=0..k)*t^k,k=0..80),t,80): # GUESSING (S2) > P:=subs(Fx0(t)=T,gfun:-seriestoalgeq(S,Fx0(t)) [START_REF] Abel | OEuvres complètes. Tome II[END_REF]): # RIGOROUS PROOF (S3) > ker := (T,t,x) -> (x+T+x^2*T^2)*t-x*T: > pol := unapply(P,T,t,x): > p1 := resultant(pol(z-T,t,x),ker(t*z,t,x),z): > p2 := subs(T=x*T,resultant(numer(pol(T/z,t,z)),ker(z,t,x),z)): > normal(primpart(p1,T)/primpart(p2,T)); 1 # time (in sec) and memory consumption (in Mb) > trunc(time()-st),trunc((kernelopts(bytesused)-bu)/1000^2); 7, 617 1. (Reduced kernel equation) Plugging

y 0 = x -t - √ -4t 2 x 3 + x 2 -2tx + t 2 2tx 2 = t + 1 x t 2 + x 3 +1 x 2 t 3 + 3x 3 +1 x 3 t 4 + 2x 6 +6x 3 +1 x 4 t 5 + • • • ∈ Q[x, x -1 ][[t]],
in [START_REF] Arquès | Dénombrements de chemins dans R 2 soumis à contraintes[END_REF] shows that U = K(t; x, 0) satisfies the reduced kernel equation

(18) 0 = x • y 0 -x • t • U(t, x) -y 0 • t • U(t, y 0 ).
2. (Uniqueness) Eq. ( 18) has a unique solution in

Q[[x, t]], namely U = K(t; x, 0). 3. (Proof) Candidate P x,0 (T, t, x) guessed in (23) admits a root H in Q[[t, x]].
Resultant computations prove that U = H(t, x) also satisfies [START_REF] Askey | How can mathematicians and mathematical historians help each other?[END_REF]. By uniqueness, K(t; x, 0) coincides with H(t, x), which is algebraic. In the case of Gessel walks the proof follows the same strategy, but several complications occur:

• the diagonal symmetry of the step set is lost, so G(t; x, y) = G(t; y, x); • G(t; 0, 0) occurs in [START_REF] Albert | Permutations sortable by two stacks in parallel and quarter plane walks[END_REF] (because of the step );

• guessed equations are ≈ 5 000 times bigger. To bypass these difficulties, one ingredient of the solution proposed in [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF] is to replace equation ( 18) by an equivalent system of two reduced kernel equations, and to make use of fast algorithms for manipulating algebraic series, inspired by the algorithms for sums and products of algebraic numbers, designed in [START_REF] Bostan | Fast computation of special resultants[END_REF]. For more details, we refer the reader to the original article [START_REF] Bostan | The complete generating function for Gessel walks is algebraic[END_REF].

D-finiteness proofs for models A and B in Fig. 13.

Here we simply state two recent results that have been discovered and proved using the guess-and-prove strategy explained before. The first one, Theorem 27, is remarkable in that it is a (more difficult) analogue of Theorem 11 (former Gessel Conjecture 1). The simple formulas beg for a combinatorial proof, but for the moment no human proof at all is known for it. Theorem 27 ([60]). The generating function E(t) = F A (t; 0, 0) = F B (t; 0, 0) of excursions for the quadrant models A and B in Fig. 13 is

4 F 3 5 6 1 2 1 7 6 5 4 2 7 4 27t 2 = ∑ n 0 6(6n + 1)!(2n + 1)! (3n)!(4n + 3)!(n + 1)! t 2n = 1 + 3t 2 + 26t 4 + 323t 6 + • • • .
It is algebraic of degree 6, root of the polynomial

16t 10 T 6 + 48t 8 T 5 + 8(6t 2 + 7)t 6 T 4 + 32(3t 2 + 1)t 4 T 3 + (48t 4 -8t 2 + 9)t 2 T 2 + (48t 4 -56t 2 + 1)T + (16t 4 + 44t 2 -1). A parametric expression of E(t) is t 2 E(t) = Z(1 -6Z + 4Z 2 )
, where Z is the unique series in t with constant term 0 satisfying

Z(1 -Z)(1 -2Z) 4 = t 2 .
The second result, Theorem 28, has two parts. The first part is remarkable in that it provides the first example of D-finiteness result of a (non-algebraic) quadrant model that is currently proved uniquely via computer algebra. The second part is remarkable in that it is a (more difficult) analogue of Theorem 12 (former Gessel Conjecture 2). Again, no human proof is known for these results. Theorem 28 ([60]). (a) The full generating function F A (t; x, y) is D-finite. (b) The full generating function F B (t; x, y) is algebraic, of degree 12. It satisfies

F B (t; x, y) = xy -t(1 + x 2 )F B (t; x, 0) -t(1 + y)F B (t; 0, y) + tF B (t; 0, 0) (y -t(1 + y)(1 + x 2 (1 + y))
) .

The sections F B (t; x, 0) and F B (t; 0, y) can be written in parametric form as follows.

Let T(t) = t + 4t 3 + 48t 5 + • • • be the unique series in t with constant term 0 such that

T(1 -4T 2 ) = t.
Let S(t) = t + 5t 3 + 62t 5 + • • • be the unique series in t with constant term 0 such that

S(1 -S 2 ) 2 = t(1 + S 2 ) 3 .
Then F B (t; x, 0) has degree 12 and is quadratic over Q(x, S):

F B (t; x, 0) = 1 + S 2 1 -S 2 3 × x(1 + 6S 2 + S 4 )-2S(1 -S 2 )(1 + x 2 )-(x -2S + xS 2 ) (1 -S 2 ) 2 -4xS(1 + S 2 ) 2x(1 + x 2 )S 2 .
Let finally W(t, y) be the unique power series in t with constant term 0 such that

W (1 -(1 + y)W) = T 2 .
Then F B (t; 0, y) has degree 6 and is rational in T and W:

F B (t; 0, y) = t -2 W(1 -4T 2 -2W).
Moreover, its coefficients are doubly hypergeometric:

F B (t; 0, y) = ∑ n j 0 6(2j + 1)!(6n + 1)!(2n + j + 1)! j! 2 (3n)!(4n + 2j + 3)!(n -j)!(n + 1)
y j t 2n .

Transcendence proofs for D-finite models.

We have seen that the guessand-prove paradigm can be successfully used to prove D-finiteness and algebraicity. The proofs are constructive by design: they internally construct (differential, or algebraic) equations. It might thus look surprising that guess-and-prove can also be used to prove transcendence, that is, lack of algebraic equations. The framework is as follows. Assume that f ∈ Q[[t]] is a D-finite power series for which some linear differential equation L( f ) = 0 (not necessarily of minimal order) is known. For instance, this differential equation could have been produced itself by a guessand-prove process. The question is how to prove that f is transcendental? This is interesting especially in cases where all known transcendence criteria (such as those in [START_REF] Flajolet | Analytic models and ambiguity of context-free languages[END_REF]) fail to apply. Such cases do occur, as seen in §1.20 for the length generating function F(t; 1, 1) for models 5-10 in Fig. 8, for which the asymptotic behavior is not incompatible with algebraicity. For these models, one possible workaround uses the factorization patterns of the differential operators for F(t; 1, 1): the operators systematically factor as a product of an order-2 operator on the left, and several order-1 operators on the right, so that Kovacic's algorithm [START_REF] Kovacic | An algorithm for solving second order linear homogeneous differential equations[END_REF] can be used to prove transcendence in an uniform way [START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF].

But factorization of linear differential operators, although quite well studied in theory [START_REF] Grigor Ev | Complexity of factoring and calculating the GCD of linear ordinary differential operators[END_REF][START_REF] Singer | Testing reducibility of linear differential operators: a group-theoretic perspective[END_REF][START_REF] Bronstein | An introduction to pseudo-linear algebra[END_REF][START_REF] Van Hoeij | Factorization of differential operators with rational functions coefficients[END_REF] is computationally very expensive, or even infeasible in practice, when applied to operators of high orders. Such an example is provided by Model B in Fig. 13. By Theorem 28, its full generating function F A (t; x, y) is D-finite, and by Theorem 27 its excursions generating function F A (t; 0, 0) is even algebraic. A natural question is: is F A (t; x, y) algebraic, or transcendental? The answer is contained in the Theorem 29 below.

Theorem 29 ([60]). F

A (t; 1, 0) = 1 + t + 4 t 2 + 8 t 3 + 39 t 4 + 98 t 5 + 520 t 6 + • • • is transcendental.
In particular, the full generating function F A (t; x, y) is transcendental.

The only available proof [START_REF] Bostan | On 3-dimensional lattice walks confined to the positive octant[END_REF] uses the guess-and-prove method. It consists in computing the minimal-order operator L f min for f = F A (t; 1, 0) and checking that L f min admits logarithms in some local expansions, which in particular prevents algebraicity of f . The computation of L f min is inspired by [380, §9]. The main idea can be traced back at least to [START_REF] Ohtsuki | On the number of apparent singularities of a linear differential equation[END_REF]; similar arguments are used in [START_REF] Chudnovsky | Rational and Padé approximations to solutions of linear differential equations and the monodromy theory[END_REF][START_REF] Bertrand | Équations différentielles linéaires et majorations de multiplicités[END_REF] and [149, §2].

All in all, the argument may be viewed as a general technique for proving transcendence of D-finite power series ; it reduces the transcendence question to dif-There exists an alternative algorithmic procedure based on [START_REF] Singer | Algebraic solutions of nth order linear differential equations[END_REF], that allows in principle to answer this question [351]. It involves, among other things, factoring linear differential operators, and deciding whether a linear differential operator admits a basis of algebraic solutions. However, this procedure would have a very high computational cost when applied to our situation. ferential guessing. In the case of f = F A (t; 1, 0), the proof consists in the following steps:

1

Definition.

Let K be a field, typically Q or a finite field F p for a prime p. Given a column vector of power series

F = ( f 1 , . . . , f n ) T ∈ K[[t]
] n and an n-tuple of integers d = (d 1 , . . . , d n ) ∈ N n , a Hermite-Padé approximant of type d for F is a row vector of polynomials P = (P 1 , . . . , P n ) ∈ K[t] n \ {0} such that:

(1)

P • F = P 1 f 1 + • • • + P n f n = O(t σ ) with σ = ∑ i (d i + 1) -1, (2)
deg(P i ) ≤ d i for all i. The integer σ is called the order of the approximant P, and d is called its type.

When n = 2, Hermite-Padé approximants are called Padé approximants, a notion intimately related to rational approximations and continued fractions. When f = A ( -1) , resp. f = A -1 , for some A ∈ K[[t]], we talk about differential approximation, resp. of algebraic approximation, which form the basis of the differential, resp. algebraic, guessing described in §2.3.2.

These concepts were initially studied by Hermite [START_REF] Hermite | Sur la généralisation des fractions continues algébriques[END_REF] and by Padé [START_REF] Padé | Sur la généralisation des fractions continues algébriques[END_REF], and turned out to be very useful in irrationality and transcendence proofs. For instance they (or, variants of them) served to prove the transcendence of e [START_REF] Hermite | Sur la fonction exponentielle[END_REF] and of π [START_REF] Lindemann | Ueber die Zahl π[END_REF], see also [START_REF] Mahler | Applications of some formulae by Hermite to the approximation of exponentials and logarithms[END_REF][START_REF] Mahler | Perfect systems[END_REF][START_REF] Beukers | Padé-approximations in number theory[END_REF]. The Chudnovsky brothers [START_REF] Chudnovsky | Rational and Padé approximations to solutions of linear differential equations and the monodromy theory[END_REF][START_REF] Chudnovsky | Padé approximation and the Riemann monodromy problem[END_REF][START_REF] Chudnovsky | Applications of Padé approximations to Diophantine inequalities in values of G-functions[END_REF] used Hermite-Padé approximants for irrationality and transcendence proofs for values of quite general D-finite functions. A spectacular recent success using such approximants is the proof [START_REF] Ball | Irrationalité d'une infinité de valeurs de la fonction zêta aux entiers impairs[END_REF] of the irrationality of infinitely many values of the zeta function at odd integers. In most of these works, arithmetic results are obtained using explicit closed-form expressions for approximants, highly based on the structure of the functions to be approximated.

Our need is different, of algorithmic nature: we need fast algorithms that compute Hermite-Padé approximants on generic inputs. Before showing how to do that, we start with a very basic example.

Worked example. Let us compute a Hermite-Padé approximant of type

(1, 1, 1) for (1, C, C 2 ), where C(t) = 1 + t + 2t 2 + 5t 3 + 14t 4 + 42t 5 + O(t 6 ). * * This boils down to finding α 0 , α 1 , β 0 , β 1 , γ 0 , γ 1 ∈ Q (not all zero) such that α 0 +α 1 t+(β 0 + β 1 t)(1 + t + 2t 2 + 5t 3 + 14t 4 )+(γ 0 + γ 1 t)(1 + 2t + 5t 2 + 14t 3 + 42t 4 ) = O(t 5
])[1][1]; / 2 \ /d \ 2 |d | 2 y(t) + (10 t -2) |--y(t)| + (4 t -t) |---y(t)| \dt / | 2 | \dt / Figure 16.
A toy guessing example: the generating function of the Catalan sequence is recognized to be algebraic by gfun using its first 5 terms. It is also recognized to be D-finite using a few more terms.

Identifying coefficients, this is equivalent to a homogeneous linear system:

      1 0 1 0 1 0 0 1 1 1 2 1 0 0 2 1 5
2 0 0 5 2 14 5 0 0 14 5 42 14

      ×         α 0 α 1 β 0 β 1 γ 0 γ 1         = 0 ⇐⇒       1 0 1 0 1 0 1 1 1 2 0 0 2 1 5 0 0 5 2 14 0 0 14 5 42       ×       α 0 α 1 β 0 β 1 γ 0       = -γ 1       0 1 2 5 14      
. By homogeneity, one can choose γ 1 = 1. The bottom-right 3 × 3 minor shows that one can take (β 0 , β 1 , γ 0 ) = (-1, 0, 0). Finally, the other values are α 0 = 1, α 1 = 0.

Thus the searched approximant is (1, -1, t): this means that we have guessed the candidate P = 1y + ty 2 such that P(t, C(t)) = 0 mod t 5 . This kind of functionality is implemented in most Computer Algebra systems. For instance, Maple's package gfun [START_REF] Salvy | Gfun: a Maple package for the manipulation of generating and holonomic functions in one variable[END_REF] implements the commands seriestoalgeq and listtoalgeq for algebraic approximants, resp. seriestodiffeq and listtodiffeq for differential approximants.

Existence and quasi-optimal computation.

The existence of Hermite-Padé approximants is guaranteed by a simple linear algebra argument: the undetermined coefficients of a potential approximant P = (∑ d i j=0 p i,j t j ) i ∈ K[t] n satisfy a linear homogeneous system with σ = ∑ i (d i + 1) -1 equations and σ + 1 unknowns. This proof is constructive and gives a first, naive, algorithm for the effective computation of Hermite-Padé approximants, of complexity O(σ ω ), where 2 ≤ ω ≤ 3 denotes a feasible linear algebra exponent, that is a constant that governs the complexity of most operations on dense matrices with coefficients in K [START_REF] Zur Gathen | Modern computer algebra[END_REF]Ch. 12].

However, as can be seen on the example in §2.5.2, the linear system has a Toeplitz-like structure: its matrix is obtained by concatenation of Sylvester-like blocks, that possess the Toeplitz property of invariance along diagonals. There are better algorithms that are able to exploit this structure. For instance, a generalization of the Euclidean algorithm yields a fast algorithm, of quadratic complexity O(σ 2 ) [START_REF] Sergeyev | A recursive algorithm for Padé-Hermite approximations[END_REF] with respect to the order of approximation σ, see also [START_REF] Beckermann | A uniform approach for the fast computation of matrix-type Padé approximants[END_REF] and the references therein. There are even faster algorithms that achieve a complexity softly-linear in σ, namely O(σ log 2 σ). They are based on fast (FFT-based) polynomial multiplication [START_REF] Zur Gathen | Modern computer algebra[END_REF]Chap. 8], and they rely on a divide-and-conquer scheme. Some are direct [START_REF] Beckermann | A uniform approach for the fast computation of matrix-type Padé approximants[END_REF], other use the artillery of the theory of matrices with small displacement rank [START_REF] Pan | Structured matrices and polynomials[END_REF][START_REF] Bostan | Solving structured linear systems with large displacement rank[END_REF][START_REF] Bostan | On Matrices With Displacement Structure: Generalized Operators and Faster Algorithms[END_REF].

Here we give a rough sketch of the structure of the superfast Beckermann-Labahn algorithm [29, §4], when applied to compute a Hermite-Padé approximant of type (d, . . . , d) for F = ( f 1 , . . . , f n ) ∈ K[[t]] n . The two main ideas are: to compute a whole matrix of approximants instead of just one approximants; to use a strategy of divide-and-conquer with respect to the order of the approximant σ = n(d + 1) -1. The algorithm proceeds as follows:

1. If σ is below some chosen threshold, then use the naive algorithm 2. Else: (a) recursively compute

P 1 ∈ K[t] n×n s.t. P 1 • F = O(t σ/2 ), deg(P 1 ) ≈ d 2 (b) compute the residue R ∈ K[[t]] n×n s.t. P 1 • F = t σ/2 • (R + O(t σ/2 )) (c) recursively compute P 2 ∈ K[t] n×n s.t. P 2 • R = O(t σ/2 ), deg(P 2 ) ≈ d 2 (d) return P := P 2 • P 1 By construction, P • F = O(t σ ).
The precise choices of degrees is a delicate issue, and is one of the most difficult technical parts in the correctness proof. From the complexity point of view, up to logarithmic factors, the total cost of the whole algorithm is concentrated into the one of a product of n × n polynomial matrices of degree ≈ d 2 , that is Õ(n ω d) operations in K. For more details, the reader is referred to the original article [START_REF] Beckermann | A uniform approach for the fast computation of matrix-type Padé approximants[END_REF].

2.6. Back to the exercise in §1.1. To finish this section, we come back to the problem stated at the very beginning of the memoir, and show how to apply the Hermite-Padé approximation in order to guess the answer. A rigorous proof will be given in §3. [START_REF] Aebly | Démonstration du problème du scrutin par des considérations géométriques[END_REF]. In what follows, S denotes the step set {↑, ←, }.

A recurrence relation for S-walks in Z × N.

Let us denote by h n;i,j the number of S-walks in Z × N of length n from (0, 0) to (i, j).

The numbers h n;i,j satisfy the following recurrence:

h n;i,j =        0 if j < 0 or n < 0, 1 i=j=0 if n = 0, ∑ (k, )∈S h n-1;i-k,j-otherwise.
The following Maple lines compute the first terms of the generating function A of the sequence (a n ) n = (h n;0,0 ) n counting S-walks in Z × N that end at the origin:

> h:=proc(n,i,j) option remember; if j<0 or n<0 then 0 elif n=0 then if i=0 and j=0 then 1 else 0 fi else h(n-1,i,j-1) + h(n-1,i+1,j) + h(n-1,i-1,j+1) fi end: > A:=series(add(h(n,0,0)*t^n, n=0..30), t,30); They produce the output A =1 +

A recurrence relation for S-walks in N

2 . Let us denote by q n;i,j the number of S-walks in N 2 of length n from (0, 0) to (i, j).

The numbers q n;i,j satisfy the same recurrence as h n;i,j , but with different boundary conditions:

q n;i,j =        0 if i < 0 or j < 0 or n < 0, 1 i=j=0 if n = 0, ∑ (k, )∈S q n-1;i-k,j-otherwise.
The following Maple lines compute the first terms of the generating function B of the sequence (b n ) n = (∑ k q n;k,k ) n , counting S-walks in N 2 that end on the diagonal:

> q:=proc(n,i,j) option remember; if i<0 or j<0 or n<0 then 0 elif n=0 then if i=0 and j=0 then 1 else 0 fi else q(n-1,i,j-1) + q(n-1,i+1,j) + q(n-1,i-1,j+1) fi end: > B:=series(add(add(q(n,k,k), k=0..n)*t^n, n=0..30), t,30);

The produced output is

B =1 + 3 t 3 +
We observe that A = B mod t 30 , but of course this is not yet a proof that A = B.

Guessing a closed form for the answer.

From the first 30 terms of A and B, one can guess a nice formula for them. The following Maple lines show a way to do that. One could first guess a differential equation (seriestodiffeq), then convert it to a recurrence (diffeqtorec); here we appeal to a shortcut (seriestorec) which guesses directly a first-order recurrence for the coefficients of A. The series A is a hypergeometric function, that can be computed explicitly.

> gfun:-seriestorec(A, u(n))[1]; 2 2 {(-27 n -81 n -54) u(n) + (n + 9 n + 18) u(n + 3), u(0) = 1, u(1) = 0, u(2) = 0} > rsolve(%, u(n)): > A:=sum(subs(n=3*n, op(2,%))*t^(3*n), n=0..infinity); 3 A := hypergeom([1/3, 2/3], [2], 27 t )
In other words, guessing predicts the following equality, equivalent to (1):

A(t) = B(t) = 2 F 1 1 3 2 3 2 27 t 3 = ∞ ∑ n=0 (3n)! n! 3 t 3n n + 1 .

Creative telescoping.

Then we wish to show that (n

+ 1) 2 b n+1 -n 2 b n-1 = (11n 2 + 11n + 3)b n , where b n = ∑ n k=0 F n,k with F n,k = ( n k ) 2 ( n+k k
). Neither Cohen nor I had been able to prove this in the intervening 2 months. After a few days of fruitless effort the specific problem was mentioned to Don Zagier (Bonn), and with irritating speed he showed that indeed the sequence (b n ) satisfies this recurrence. We cleverly construct B n,k = (k

2 + 3(2n + 1)k -11n 2 -9n -2)F n,k , with the motive that B n,k -B n,k-1 = (n + 1) 2 F n+1,k -(11n 2 + 11n + 3)F n,k -n 2 F n-1,k ,
and, O mirabile dictu, the sequence (b n ) does indeed satisfy the recurrence by virtue of the method of creative telescoping.

A. van der Poorten [START_REF] Van Der Poorten | Apéry's proof of the irrationality of ζ(3)[END_REF].

3.1. Diagonals. Algebraic power series are D-finite ( §1.11). An intermediate important class of power series is formed by diagonals of rational functions. All the examples of D-finite generating functions occurring in our combinatorial context of enumeration of walks appear to be diagonals, either directly (by their combinatorial definition), or indirectly (by the resolution method). The differential equations that they satisfy are special cases of Picard-Fuchs equations for periods of rational functions, and can be constructed algorithmically. A conjecture of Christol's [START_REF] Christol | Globally bounded solutions of differential equations[END_REF] predicts even more: any D-finite power series S ∈ Z[[t]] with finite non-zero radius of convergence is the diagonal of a rational function.

In combinatorics, the importance of diagonals stems from the fact that numerous combinatorial constructions on generating functions (Hadamard products, constant terms or positive parts of Laurent series, . . . ) can be encoded as diagonals [357]. A classical result [START_REF] Lipshitz | The diagonal of a D-finite power series is D-finite[END_REF][START_REF] Christol | Diagonales de fractions rationnelles[END_REF] asserts that the diagonal of a rational function is D-finite (Theorem 33). A natural question is then: how to obtain algorithmically the linear differential equation satisfied by a diagonal? The problem can be reformulated in terms of the computation of a multiple integral with parameters taken on a cycle ( §3), and can thus be attacked from a geometrical viewpoint [START_REF] Beukers | Irrationality of π 2 , periods of an elliptic curve and Γ 1 (5)[END_REF][START_REF] Deligne | Intégration sur un cycle évanescent[END_REF].

∈ Q[[x 1 , . . . , x n ]] F = ∑ i 1 ,...,i n ≥0 a i 1 ,...,i n x i 1 1 • • • x i n n
is the univariate power series Diag(F) = ∑ i a i,...,i t i .

Pólya's theorem.

Almost a century ago, Pólya proved that diagonals of bivariate rational functions are algebraic [START_REF] Pólya | Sur les séries entières, dont la somme est une fonction algébrique[END_REF]. Later, Furstenberg showed that the converse also holds [START_REF] Furstenberg | Algebraic functions over finite fields[END_REF]. Interestingly, Pólya's result becomes false for more than two variables. A simple example is provided by

Diag 1 1 -x -y -z = ∑ n≥0 3n n, n, n t n = 2 F 1 1 3 2 3 1 27t .
The transcendence of this series can be proved in various ways, for instance by using the asymptotics

( 3n n,n,n ) = (3n)! n! 3 ∼ 3 3n √ 3
2πn and [188, Theorem D], or by using the interlacing criterion from Theorem 8.

Pólya's result can be proved as follows. First, using the simple observation Diag(F)(t) = [x 0 ] F(x, t/x), the diagonal of the rational function F(x, y) ∈ Q(x, y) is encoded as a complex integral using Cauchy's integral theorem (for some > 0)

Diag (F) (t) = [x -1 ] 1 x F x, t x = 1 2πi |x|= F x, t x dx x ,
which in a second step can be evaluated using the residues theorem as a sum of residues (precisely: the residues of F(x, t/x)/x at its "small poles", having limit 0 at t = 0). Each of these residues are algebraic functions, and so is their sum Diag(F).

Example 31 (Dyck bridges). Let S = { , }, and let B n be the number of Dyck bridges (i.e. S-walks in Z 2 starting at (0, 0) and ending on the horizontal axis), of length n. Using a rotation counterclockwise by π/4, the integer B n is seen to be the number of {↑, →}-walks in Z 2 from (0, 0) to (n, n). This implies

B(t) = ∑ n≥0 B n t n = Diag 1 1 -x -y ,
and the proof sketched above concludes:

B(t) = 1 2πi |x|= dx x -x 2 -t = 1 1 -2x x= 1-√ 1-4t 2 = 1 √ 1 -4t = ∑ n≥0 2n n t n .
Rothstein-Trager resultant. It is not always possible to compute explicitly a closedform expression for the poles and the residues, as we did in Example 31, for instance when the denominator of F(x, t/x)/x has degree more than 4. However, using resultants one can compute annihilating polynomials for them, and thus also for the diagonal. We show how this is done if F(x, t/x)/x has simple poles only. Assume that K is a field (in our case, K is a placeholder for Q(t)), and let A, B ∈ K[x] be such that deg(A) < deg(B), with B squarefree. Then the rational function F = A/B has simple poles only, and if F admits the partial fraction decomposition

F = ∑ i γ i x -β i , then the residue γ i of F at the pole β i equals γ i = A(β i ) B (β i )
. Therefore, the residues γ i of F are roots of the so-called Rothstein-Trager resultant [START_REF] Rothstein | Aspects of symbolic integration and simplification of exponential and primitive functions[END_REF][START_REF] Trager | Algebraic factoring and rational function integration[END_REF]:

R(t) = Res x B(x), A(x) -t • B (x) ,
which was originally introduced in computer algebra for the symbolic (indefinite) integration of rational functions.

A generalization of the Rothstein-Trager resultant to the case of multiple poles was given by Bronstein [START_REF] Bronstein | Formulas for series computations[END_REF]. Example 32 (Diagonal Rook paths). Consider the following question [START_REF] Erickson | Enumerating rook and queen paths[END_REF][START_REF] Coker | Enumerating a class of lattice paths[END_REF]: A Rook can move any number of squares horizontally or vertically. Assuming that the Rook moves right or up at each step, how many paths can the Rook take from the lower-left corner square to the upper-right corner square of an N × N chessboard? Denote this number by d N , see Fig. 18. The generating function of (d n ) n is

Diag(F) = [x 0 ] F(x, t/x) = 1 2πi F(x, t/x) dx x , where F = 1 1 -x 1-x -y 1-y .
Then, Diag(F) is a sum of roots y(t) of the Rothstein-Trager resultant > F:=1/(1-x/(1-x)-y/(1-y)): > G:=normal(1/x*subs(y=t/x,F)): > factor(resultant(denom(G),numer(G)-y*diff(denom(G),x),x));

which is t 2 (1 -t)(2y -1)(36ty 2 -4y 2 + 1 -t).
By identifying which residues correspond to small poles, one concludes that the generating function of diagonal 2D

Rook paths is equal to the algebraic function

1 2 1 + 1 -t 1 -9t .
The same method can be used for other walks of the same type [START_REF] Kung | Catalan lattice paths with rook, bishop and spider steps[END_REF]. Algorithmic questions related to the computation of algebraic equations for diagonals Diag(F) of bivariate rational functions F have been considered in connection with the enumeration of 1D lattice walks (bridges, excursions and meanders) by Banderier and Flajolet [START_REF] Banderier | Basic analytic combinatorics of directed lattice paths[END_REF]. A general and efficient algorithm that computes an annihilating polynomial for Diag(F) was later proposed by Bostan, Dumont and Salvy [START_REF] Bostan | Algebraic diagonals and walks[END_REF][START_REF] Bostan | Algebraic diagonals and walks: algorithms, bounds, complexity[END_REF]; that solves positively the question of an effective version of Pólya's theorem. On the negative side, they showed that the minimal polynomial of Diag(F) has generically an exponential size with respect to the degree of the input rational function F. By contrast, linear differential equations satisfied by Diag(F) had been proved to have polynomial size [START_REF] Bostan | Complexity of creative telescoping for bivariate rational functions[END_REF]. This implies that for bivariate diagonals, the differential equations are the right data-structure, and not algebraic equations. It was shown in [START_REF] Bostan | Algebraic diagonals and walks[END_REF][START_REF] Bostan | Algebraic diagonals and walks: algorithms, bounds, complexity[END_REF] that the first N terms of generating functions for various 1D walks can be computed in quasi-linear complexity in N using this data-structure.

Lipshitz's theorem.

Although diagonals of multivariate rational functions are not necessarily algebraic, they are still D-finite. In fact, much more holds.

Theorem 33 (Lipshitz,[START_REF] Lipshitz | The diagonal of a D-finite power series is D-finite[END_REF]). Diagonals of D-finite power series are D-finite.

For rational † † ‡ ‡ power series, this result was previously obtained by Christol in an "elementary" way under a regularity assumption [START_REF] Christol | Diagonales de fractions rationnelles et équations différentielles[END_REF], and in the general case using quite involved geometric arguments [START_REF] Deligne | Intégration sur un cycle évanescent[END_REF][START_REF] Christol | Diagonales de fractions rationnelles et équations de Picard-Fuchs[END_REF][START_REF] Christol | Diagonales de fractions rationnelles[END_REF], see also [START_REF] Lipshitz | Rational functions, diagonals, automata and arithmetic[END_REF][START_REF] Christol | Globally bounded solutions of differential equations[END_REF][START_REF] Christol | Diagonals of rational fractions[END_REF]. Very briefly, the argument is the following. First, as in the bivariate case, if f ∈ Q(x 1 , . . . , x n ) ∩ Q[[x 1 , . . . , x n ]], the residue theorem allows to write (for some > 0) ( 21)

Diag( f )(t) = 1 (2πi) n-1 |x 1 |=•••=|x n-1 |= f x 1 , . . . , x n-1 , t x 1 • • • x n-1 dx 1 • • • dx n-1 x 1 • • • x n-1 ,
so that Diag( f )(t) is the (relative) period of a (family of) rational function(s) [START_REF] Katz | On the differential equations satisfied by period matrices[END_REF]. Its D-finiteness is then a consequence of the (highly non-trivial) finite-dimension property over C(t) of the de Rham cohomology for the complement of the variety in A n C(t) defined by the equations denom( f )(x 1 , . . . , x n ) = 0 and x 1 • • • x n = t. § § In more down-to-earth terms this proof guarantees, in a non-effective way, that repeated differentiation under the integral sign eventually produces a finite sequence of rational integrands that admit a linear combination with coefficients in Q(t) that becomes an exact differential. On the one hand, this geometric method allows access to more information about the minimal-order differential equation: it is Fuchsian and it has only rational exponents at each singularity (see [START_REF] Griffiths | On the periods of integrals on algebraic manifolds[END_REF][START_REF] Griffiths | On the periods of certain rational integrals. I, II[END_REF][START_REF] Griffiths | Periods of integrals on algebraic manifolds: Summary of main results and discussion of open problems[END_REF] for an analytic proof and [START_REF] Katz | Nilpotent connections and the monodromy theorem: Applications of a result of Turrittin[END_REF][START_REF] Katz | The regularity theorem in algebraic geometry[END_REF] for an arithmetic proof ¶ ¶ ). On the other hand, it is non constructive. (See §3.2.3 for a way to make it constructive, using the so-called Griffiths-Dwork reductions.) By contrast, Lipshitz' proof [START_REF] Lipshitz | The diagonal of a D-finite power series is D-finite[END_REF] is elementary and constructive. However, the algorithm behind its proof in highly inefficient. We demonstrate this using the example provided by the following combinatorial problem.

Example 34 (Diagonal 3D Rook paths). Consider the following question [START_REF] Erickson | Enumerating rook and queen paths[END_REF]: How many ways can a Rook move on a 3D chessboard from (0, 0, 0) to (N, N, N), where each step is a positive integer multiple of (1, 0, 0), (0, 1, 0), or (0, 0, 1)? This is a 3D extension of the question in Example 32. Denote by D N the number of diagonal 3D Rook paths of length N. The first terms of the sequence (D N ) are:

1, 6, 222, 9918, 486924, 25267236, 1359631776, 75059524392, . . . The combinatorial problem readily translates into an algebraic problem. The generating function ∆(t) = ∑ n≥0 D n t n of diagonal 3D Rook paths is the diagonal of † † For rational series, Th. [START_REF] Lipshitz | The diagonal of a D-finite power series is D-finite[END_REF] had been conjectured by Stanley [356,§4(b)] and incompletely proved in [START_REF] Zeilberger | Sister Celine's technique and its generalizations[END_REF][START_REF] Gessel | Two theorems on rational power series[END_REF].

‡ ‡ For algebraic series, Th. [START_REF] Lipshitz | The diagonal of a D-finite power series is D-finite[END_REF] can be proved by reduction to the rational case [START_REF] Denef | Algebraic power series and diagonals[END_REF][START_REF] Adamczewski | Diagonalization and rationalization of algebraic Laurent series[END_REF], for the price of doubling the number of variables.

§ § The finiteness proof needs Hironaka's resolution of singularities, among other things [START_REF] Grothendieck | On the de Rham cohomology of algebraic varieties[END_REF][START_REF] Monsky | Finiteness of de Rham cohomology[END_REF][START_REF] Hartshorne | On the de Rham cohomology of algebraic varieties[END_REF]]. ¶ ¶ Katz first shows in [246, §5] that the minimal-order equation for a period is globally nilpotent; this result has been generalized by the Chudnovskys to any G-function [START_REF] Chudnovsky | Applications of Padé approximations to Diophantine inequalities in values of G-functions[END_REF], see also [START_REF] Dwork | An introduction to G-functions[END_REF]Chap. VIII]. Then, Katz shows in [246, §13] that globally nilpotent operators are Fuchsian, with rational exponents; see also [START_REF] Honda | Algebraic differential equations[END_REF][START_REF] Dwork | Differential operators with nilpotent p-curvature[END_REF] for a more elementary proof. the rational function F(x, y, z) given by

1 -∑ n≥1 x n -∑ n≥1 y n -∑ n≥1 z n -1 = (1 -x)(1 -y)(1 -z) 1-2(x+y+z)+3(xy+yz+zx)-4xyz .
A closed form for ∆(t) has been obtained by Bostan, Chyzak, van Hoeij and Pech [START_REF] Bostan | Explicit formula for the generating series of diagonal 3D rook paths[END_REF], as an integral of a hypergeometric 2 F 1 . Its form is very similar to the ones in §1.19 for quarter-plane walks with small steps (models 1-19 in Fig. 7).

Theorem 35 ([77]).

∑ n D n t n = 1 + 6 • t 0 2 F 1 1/3 2/3 2 27x(2-3x) (1-4x) 3 (1 -4x)(1 -64x) dx.
The proof of Theorem 35 consists in first computing a differential equation for Diag(F), then in solving it in closed form, using algorithms in [START_REF] Bostan | Explicit formula for the generating series of diagonal 3D rook paths[END_REF][START_REF] Kunwar | Second order differential equations with hypergeometric solutions of degree three[END_REF][START_REF] Imamoglu | Computing hypergeometric solutions of second order linear differential equations using quotients of formal solutions[END_REF][START_REF] Imamoglu | Computing hypergeometric solutions of second order linear differential equations using quotients of formal solutions and integral bases[END_REF]. In what follows, we focus on the first part, and describe the main steps on Lipshitz' proof when applied to prove the D-finiteness of Diag(F). The starting point is the following: If one is able to find a nonzero differential operator of the form

L(t, ∂ t , ∂ x , ∂ y ) = P(t, ∂ t ) + ( higher-order terms in ∂ x and ∂ y ) that annihilates G = 1 xy • F x, y x , t y
, then P(t, ∂ t ) annihilates Diag(F). This is explained by the sequence of equalities:

0 = [x -1 y -1 ]L(G) = [x -1 y -1 ]P(G) = P([x -1 y -1 ]G) = P(Diag(F)).
The first equality comes from 0 = L(G), the second one is a consequence of L(G) = P(G) + ∂ x (•) + ∂ y (•) and of the fact that derivatives w.r.t. x (resp. y) do not contain terms in x -1 (resp. in y -1 ); the third equality is explained by the fact that P does only depend on t; the last one comes from Diag(F) = [x 0 y 0 ] F (x, y/x, t/y). The remaining task is to show that such an L does indeed exist. To do this, a combinatorial argument is applied: By Leibniz's rule, the ( N+4 4 ) rational functions

t i ∂ j t ∂ k x ∂ y (G), 0 i + j + k + ≤ N
are contained in the Q-vector space of dimension ≤ 18(N + 1) 3 spanned by

t i x j y k denom(G) N+1 , 0 i 2N + 1, 0 j 3N + 2, 0 k 3N + 2.
Thus, if ( N+4 4 ) > 18(N + 1) 3 , then there exists L(t, ∂ t , ∂ x , ∂ y ) (resp. P(t, ∂ t )) of total degree at most N, such that LG = 0 (resp. P(Diag(F)) = 0).

At this point, note that N = 425 is the smallest integer satisfying ( N+4 4 ) > 18(N + 1) 3 . Therefore, finding the operator P by Lipshitz' argument would require solving a linear system with 1,391,641,251 unknowns and 1,391,557,968 equations! * * *

The conclusion is that Lipshitz's approach is not sufficient to obtain effectively differential equations for diagonals. This lack of efficiency motivates the creative telescoping algorithms described in the next section §3.

Creative telescoping for sums and integrals.

Toutes les relations mentionnées ci-dessus, y compris l'extraordinaire récurrence d'Apéry, sont retrouvées de manière systématique et automatique, et l'on dispose d'un outil qui permet de découvrir et de démontrer des identités d'un certain type. Le jour est sans doute proche où les formulaires classiques sur les fonctions spéciales seront remplacés par un logiciel d'interrogation performant.

P. Cartier [START_REF] Cartier | Démonstration "automatique" d'identités et fonctions hypergéométriques (d'après D. Zeilberger)[END_REF]. Creative telescoping is an algorithmic paradigm for proving identities on multiple definite integrals and sums with parameters. This powerful computer algebra tool was introduced in the early 1990s by Zeilberger in the hypergeometric/hyperexponential setting [START_REF] Zeilberger | A fast algorithm for proving terminating hypergeometric identities[END_REF][START_REF] Zeilberger | A holonomic systems approach to special functions identities[END_REF][START_REF] Almkvist | The method of differentiating under the integral sign[END_REF][START_REF] Zeilberger | The method of creative telescoping[END_REF][START_REF] Wilf | An algorithmic proof theory for hypergeometric (ordinary and "q") multisum/integral identities[END_REF], vastly generalized by Chyzak in the 2000s to the framework of holonomic functions [START_REF] Chyzak | Fonctions holonomes en calcul formel[END_REF][START_REF] Chyzak | Gröbner bases, symbolic summation and symbolic integration[END_REF][START_REF] Chyzak | Non-commutative elimination in Ore algebras proves multivariate identities[END_REF][START_REF] Chyzak | An extension of Zeilberger's fast algorithm to general holonomic functions[END_REF], and greatly enhanced and used in computerized proofs of difficult combinatorial applications by Koutschan in the 2010 [START_REF] Koutschan | A fast approach to creative telescoping[END_REF][START_REF] Koutschan | Eliminating human insight: an algorithmic proof of Stembridge's TSPP theorem[END_REF][START_REF] Koutschan | Proof of George Andrews's and David Robbins's q-TSPP conjecture[END_REF][START_REF] Koutschan | Creative telescoping for holonomic functions[END_REF][START_REF] Koutschan | Lattice Green functions of the higher-dimensional face-centered cubic lattices[END_REF][START_REF] Hassani | Lattice Green functions: the d-dimensional face-centered cubic lattice[END_REF]. Since its birth, almost 30 years ago, the methodology of creative telescoping gained more and more popularity. As of 2017, it is one of the main topics in influential conferences like ISSAC † † † , where it has yearly its own dedicated special session.

We will give a brief account on creative telescoping, since several excellent texts already exist on this topic. We refer the reader to Chyzak's habilitation thesis [START_REF] Chyzak | The ABC of Creative Telescoping -Algorithms[END_REF], and to the surveys [START_REF] Koutschan | Creative telescoping for holonomic functions[END_REF][START_REF] Chen | Some open problems related to creative telescoping[END_REF].

Example 36. (Hypergeometric summation) Creative telescoping can automatically prove the following identities: [START_REF] Dixon | On the sum of the cubes of the coefficients in a certain expansion by the binomial theorem[END_REF][START_REF] Askey | How can mathematicians and mathematical historians help each other?[END_REF])

• ∑ k∈Z (-1) k a + b a + k a + c c + k b + c b + k = (a + b + c)! a!b!c! (Dixon 1891
• a n = n ∑ k=0 n k 2 n + k k 2
satisfies the recurrence (Apéry [START_REF] Apéry | Irrationalité de ζ(2) et de ζ(3)[END_REF][START_REF] Van Der Poorten | Apéry's proof of the irrationality of ζ(3)[END_REF]) [START_REF] Strehl | Recurrences and Legendre transform[END_REF][START_REF] Strehl | Recurrences and Legendre transform[END_REF][START_REF] Strehl | Binomial identities-combinatorial and algorithmic aspects[END_REF][START_REF] Schmidt | Legendre transforms and Apéry's sequences[END_REF])

(n + 1) 3 a n+1 = (2n + 1)(17n 2 + 17n + 5)a n -n 3 a n-1 • n ∑ k=0 n k 2 n + k k 2 = n ∑ k=0 n k n + k k k ∑ j=0 k j 3 (Strehl
Example 37. (Diagonals and integrals) Creative telescoping can automatically prove the following integral and diagonal evaluations: [START_REF] Apéry | Sur certaines séries entières arithmétiques[END_REF][START_REF] Straub | Multivariate Apéry numbers and supercongruences of rational functions[END_REF])

• Diag 1 (1 -x -y)(1 -z -u) -xyzu = ∑ n 0 a n t n (Straub
• 1 2πi (1 + 2xy + 4y 2 ) exp 4x 2 y 2 1+4y 2 y n+1 (1 + 4y 2 ) 3 2 dy = H n (x) n/2 ! (Doetsch [159]) • +∞ 0 xJ 1 (ax)I 1 (ax)Y 0 (x)K 0 (x) dx = - ln(1 -a 4 ) 2πa 2 (Glasser-Montaldi [206])
where J 1 , Y 0 are Bessel functions, I 1 , K 0 are modified Bessel functions [2, Chap. 9], and H n are Hermite polynomials [START_REF] Abramowitz | Handbook of mathematical functions with formulas, graphs, and mathematical tables[END_REF]Chap. 22].

We briefly discuss the main principles of the Creative Telescoping paradigm for sums and integrals. We want to prove that I n = 2 n . The idea is that if one writes Pascal's triangle identity under the "telescopic form":

n + 1 k = n k + n k -1 = 2 n k + n k -1 - n k ,
then summation over k yields the recurrence

I n+1 = 2I n .
Taking into account the initial condition I 0 = 1 concludes the proof that I n = 2 n . More generally, assume that (u n,k ) is a bivariate sequence, and that one wants to "compute" its definite sum

F n = ∑ k u n,k ,
where "computing F n " means, as in the example, finding a recurrence relation on it. The principle is the same as in the example. Let us denote by S n and S k the forward shift operators with respect to n and k, which act on bivariate sequences by the simple rules

S n • v n,k = v n+1,k , S k • v n,k = v n,k+1
. Then, if one knows recurrence operator P(n, S n ) free of k and another recurrence operator Q(n, k, S n , S k ) such that [START_REF] Bailey | Experimental mathematics in action[END_REF] (P(n,

S n ) -(S k -1)Q(n, k, S n , S k )) • u n,k = 0,
then the sum "telescopes", leading (under "nice" boundary assumptions) to the recurrence P(n, S n ) • F n = 0. Observe that essentially the same idea was used in Lipshitz' proof of Theorem 33. The operator P is called a telescoper for u n,k , while the operator Q is called a certificate.

The whole game is then to be able to produce an equality like [START_REF] Bailey | Experimental mathematics in action[END_REF]. This is the objective of creative telescoping, a name seemingly coined by A. van der Poorten in his account of Apery's proof of the irrationality of ζ(3) [START_REF] Van Der Poorten | Apéry's proof of the irrationality of ζ(3)[END_REF], where Zagier is credited for having solved [START_REF] Bailey | Experimental mathematics in action[END_REF] for the sequence

u n,k = ( n k ) 2 ( n+k k ) 2 .
A decade later, it was Zeilberger who systematized, generalized, and quantified "Zagier's trick" in a series of articles [START_REF] Zeilberger | Sister Celine's technique and its generalizations[END_REF][START_REF] Zeilberger | A holonomic systems approach to special functions identities[END_REF][START_REF] Zeilberger | A fast algorithm for proving terminating hypergeometric identities[END_REF][START_REF] Zeilberger | The method of creative telescoping[END_REF][START_REF] Wilf | An algorithmic proof theory for hypergeometric (ordinary and "q") multisum/integral identities[END_REF]. The article [START_REF] Nemes | How to do Monthly problems with your computer[END_REF] and the entire book [325] are devoted to popularize this summation framework.

Building on previous work by Fasenmyer [START_REF] Fasenmyer | Some generalized hypergeometric polynomials[END_REF] and Verbaeten [START_REF] Verbaeten | The automatic construction of pure recurrence relations[END_REF], Wilf and Zeilberger [START_REF] Wilf | An algorithmic proof theory for hypergeometric (ordinary and "q") multisum/integral identities[END_REF] proved that the existence of a non-trivial solution (P, Q) of ( 22) is guaranteed if the summand sequence (u n,k ) is proper hypergeometric, i.e., of the form:

u n,k = p(n, k)α n β k L ∏ =1 (a n + b k + c )! ±1 , where p(n, k) ∈ Q[n, k],
where a , b , c ∈ Z, and α, β ∈ Q. Moreover, they described an algorithm to compute such a pair (P, Q), similar in spirit to that of Theorem 33, and they extended these results to multiple sums and integrals. Although based on linear algebra only, the resulting algorithm suffered from too high a complexity and from long running times in implementations, just as in the case of Lipshitz's approach for diagonals ( §3.1.2).

In parallel, Zeilberger came up with a fast algorithm for definite hypergeometric summation [START_REF] Zeilberger | A fast algorithm for proving terminating hypergeometric identities[END_REF][START_REF] Zeilberger | The method of creative telescoping[END_REF], which is based on Gosper's decision algorithm for the indefinite summation of hypergeometric sequences [START_REF] Gosper | Decision procedure for indefinite hypergeometric summation[END_REF]. Zeilberger actually realized that if the telescoper P in [START_REF] Bailey | Experimental mathematics in action[END_REF] were known beforehand, then the sequence

v n,k = Q(n, k, S n , S k ) • u n,k , which satisfies P(n, S n ) • u n,k = v n,k+1 -v n,k
, could be obtained by simply calling Gosper's algorithm. To turn this remark into an algorithm, he explained that the simultaneous search for the coefficients of the telescoper P(n, S n ) and for the rational function v n,k /u n,k amounts to using a parametrized variant of Gosper's algorithm. Zeilberger named his fast algorithm the method of creative telescoping. It is implemented in many Computer Algebra systems. In Maple, a summation package SumTools[Hypergeometric] contains a command called Zeilberger.

Example 38 (Back to the SIAM flea). Keeping notation from §1.7, the probability p n (ε) of occupying the origin at step 2n is equal to p n (ε) = ∑ n k=0 U n,k (ε), where

U n,k (ε) := 2n 2k 2k k 2n -2k n -k 1 4 + ε k 1 4 -ε k 1 4 2n-2k .
A linear recurrence for (p n (ε)) n can then be computed using Zeilberger's algorithm whose output is [START_REF] Bailey | Ten problems in experimental mathematics[END_REF] 4(n + 2) 2 S 2 n + (2n + 3) 2 (8ε 2 -1)S n + 16ε 4 (2n + 3)(2n + 1).

The probability p(ε) is equal to 1 -1 R ε (1) , where R ε (t) = ∑ n p n (ε)t n . Now (23) converts into a second-order differential equation satisfied by R ε (t), which is solved in terms of 2 F 1 's, giving the answer announced in §1.7.

Creative Telescoping for integrals.

A similar discussion applies to the case of parametrized integrals. Assume that H(t, x) is a bivariate function, and that one wants to "compute" its definite integral

I(t) = γ H(t, x) dx,
where "computing I(t)" means finding a linear differential equation satisfied by I(t). The principle from the discrete case applies to the continuous analogue. Let us denote by ∂ t and ∂ x the operators of partial derivation with respect to t and x, which act on bivariate functions by the simple rules

∂ t • f (x, t) = ∂ f ∂t and ∂ x • f (x, t) = ∂ f ∂x .
Then, if one knows a differential operator P(t, ∂ t ) free of x and another differential operator

Q(t, x, ∂ t , ∂ x ) such that (24) (P(t, ∂ t ) -∂ x Q(t, x, ∂ t , ∂ x )) • H(t, x) = 0,
then the integral with respect to x "telescopes", leading (under "nice" assumptions on the integration domain) to the differential equation

P(t, ∂ t ) • I(t) = 0.
Again, the differential operator P is called a telescoper for the integrand H(t, x), while the operator Q is called a certificate. Again, the whole game is then to be able to produce an equality like [START_REF] Bailey | Experimental mathematics in action[END_REF]. First, the existence of a pair (P, Q) like in ( 24) is guaranteed if the integrand H(t, x) is hyperexponential, that is such that both ∂H ∂t and ∂H ∂x are rational functions in t and x [START_REF] Wilf | An algorithmic proof theory for hypergeometric (ordinary and "q") multisum/integral identities[END_REF]. Moreover, the computation of such a pair (P, Q) can be done in a slow fashion, à la Lipshitz, but also by an analogue of Zeilberger's fast creative telescoping, due to Alkmkvist and Zeilberger [START_REF] Almkvist | The method of differentiating under the integral sign[END_REF]. The algorithm from [START_REF] Almkvist | The method of differentiating under the integral sign[END_REF] is implemented for instance in Maple in the DEtools package, under the name Zeilberger.

Example 39 (Diagonal Rook paths, cont.). Using notation from Example 32, one needs to compute

Diag(F) = [x 0 ] F(x, t/x) = 1 2πi F(x, t/x) dx x , where F = 1 1 -x 1-x -y 1-y .
A linear differential equation for Diag(F) can be computed using creative telescoping > F:=1/(1-x/(1-x)-y/(1-y)): > G:=normal(1/x*subs(y=t/x,F)): > DEtools [Zeilberger](G, t, x, Dt) [START_REF] Abel | OEuvres complètes. Tome II[END_REF];

whose output is (9t 2 -10t + 1)∂ 2 t + (18t -14)∂ t
and which can be solved explicitly, giving the answer Diag(F) = 1 2 1 + 1t 1 -9t .

Principle of Creative

Telescoping for multiple integrals.. In the multivariate case, we restrict our attention to the integration of rational functions. This will be sufficient for our purposes in the combinatorial applications to lattice path enumeration. Let H(t, x) be a rational function, where x = x 1 , . . . , x n denote the integration variables, and t is the parameter left after integration. Let γ be an integration domain in C n , without boundary (more precisely, an n-cycle), on which H is assumed to take finite values only. The aim is to "compute" the parametrized integral, called period,

I(t) = γ H(t, x)dx.
Example 40. The generating function for the Apéry numbers (sequence (a n ) n in Example 36) is the period of the rational integral [START_REF] Beukers | Irrationality of π 2 , periods of an elliptic curve and Γ 1 (5)[END_REF][START_REF] Beukers | A family of K3 surfaces and ζ(3)[END_REF]

1 (2πi) 3 γ dx dy dz 1 -(1 -xy)z -txyz(1 -x)(1 -y)(1 -z) ,
where γ is a suitable 3-cycle in C 3 .

More generally, diagonals of rational functions are periods, due to equation [START_REF] Bailey | Exploratory experimentation and computation[END_REF].

It is a classical theorem that periods of rational integrals are D-finite; this generalizes Theorem 33. The corresponding linear differential equations are known under the name of Picard-Fuchs equations. They describe the variation of the family of periods with respect to the parameter of the family.

Particular cases of this theorem have been proved by Legendre (1825) and Kummer (1836) [272, §29], see also [START_REF] Gray | Fuchs and the theory of differential equations[END_REF], for the periods of the complete elliptic integrals, and by Fuchs [START_REF] Fuchs | Die Periodicitätsmoduln der hyperelliptischen Integrale als Functionen eines Parameters aufgefasst[END_REF] and Picard [START_REF] Picard | Théorie des fonctions algébriques de deux variables indépendantes[END_REF] for the periods of hyperelliptic integrals and other abelian integrals on curves of arbitrary genus. The more general statements are due to Manin [START_REF] Manin | Algebraic curves over fields with differentiation[END_REF][START_REF] Manin | Rational points of algebraic curves over function fields[END_REF], who coined the term Picard-Fuchs equations (see also [START_REF] Katz | On the differential equations satisfied by period matrices[END_REF][START_REF] Katz | On the differentiation of de Rham cohomology classes with respect to parameters[END_REF][START_REF] Katz | Nilpotent connections and the monodromy theorem: Applications of a result of Turrittin[END_REF]). and to Griffiths [START_REF] Griffiths | On the periods of integrals on algebraic manifolds[END_REF][START_REF] Griffiths | On the periods of certain rational integrals. I, II[END_REF][START_REF] Griffiths | Periods of integrals on algebraic manifolds: Summary of main results and discussion of open problems[END_REF]. Modern proofs of this Dfiniteness result are based, as in the case of diagonals ( §3.1) on the finiteness of the relative de Rham cohomology of the complementary of the hypersurface defined by the singular locus of the rational integrand [START_REF] Grothendieck | On the de Rham cohomology of algebraic varieties[END_REF][START_REF] Monsky | Finiteness of de Rham cohomology[END_REF][START_REF] Hartshorne | On the de Rham cohomology of algebraic varieties[END_REF].

The principle of creative telescoping for the computation of Picard-Fuchs equations, already used by Manin in [START_REF] Manin | Rational points of algebraic curves over function fields[END_REF], is the following: if one knows a differential operator P(t, ∂ t ) free of x and rational functions (A 1 , . . . , A n ) such that ( 25)

P(t, ∂ t ) • H(t, x) = n ∑ i=1 ∂A i ∂x i ,
then the integral with respect to x "telescopes", leading to the differential equation

P(t, ∂ t ) • I(t) = 0.
(The reason is simply that integrals over cycles of pure derivatives are equal to zero.)

The differential operator P is called a telescoper for the integrand H(t, x), and (A 1 , . . . , A n ) is called a certificate. The question is then how to produce effectively an equality like [START_REF] Banderier | Limit laws for basic parameters of lattice paths with unbounded jumps[END_REF]. Ideally, one would like to compute the telescoper without computing the certificate, for reasons that will become apparent in the next example.

Example 41 (Perimeter of an ellipse). Computations of differential equations for periods can be traced back to Euler [175, §7], in his study of the perimeter p(e) of an ellipse with semi-major axis 1, as a function of its eccentricity e: The question can be casted into the framework of periods of rational integrals:

p(e) = 4 1 0 1 -e 2 x 2 1 -x 2 dx = 2π - π 2 e 2 - 3π 32
p(e) = dxdy 1 -1-e 2 x 2 (1-x 2 )y 2
, and a telescopic relation of type [START_REF] Banderier | Limit laws for basic parameters of lattice paths with unbounded jumps[END_REF] reads:

(e -e 3 )∂ 2 e + (1 -e 2 )∂ e + e •   1 1 -1-e 2 x 2 (1-x 2 )y 2   = ∂ x - e(-1-x+x 2 +x 3 )y 2 (-3+2x+y 2 +x 2 (-2+3e 2 -y 2 )) (-1+y 2 +x 2 (e 2 -y 2 )) 2 + ∂ y 2e(-1+e 2 )x(1+x 3 )y 3 (-1+y 2 +x 2 (e 2 -y 2 )) 2 .
From there, Euler's equation (ee 3 )p (e) + (1e 2 )p (e) + ep(e) = 0 follows directly. The size of the certificate is much bigger than the size of the telescoper.

Several generations of Creative Telescoping algorithms. Algorithms for creative telescoping for periods can be divided into four generations. Algorithms from the first generation (1G) -à la Lipshitzuse holonomy theory and elimination for operator ideals [START_REF] Zeilberger | A holonomic systems approach to special functions identities[END_REF][START_REF] Takayama | GrÖbner basis, integration and transcendental functions[END_REF][START_REF] Wilf | An algorithmic proof theory for hypergeometric (ordinary and "q") multisum/integral identities[END_REF][START_REF] Takayama | An approach to the zero recognition problem by Buchberger algorithm[END_REF][START_REF] Chyzak | Non-commutative elimination in Ore algebras proves multivariate identities[END_REF]; they are not very efficient in practice. Algorithms from the second generation (2G), due to Chyzak [START_REF] Chyzak | An extension of Zeilberger's fast algorithm to general holonomic functions[END_REF] and to Koutschan [START_REF] Koutschan | A fast approach to creative telescoping[END_REF], are generalizations of Zeilberger's fast algorithms for hypergeometric summation and hyperexponential integration [START_REF] Zeilberger | A fast algorithm for proving terminating hypergeometric identities[END_REF][START_REF] Zeilberger | The method of creative telescoping[END_REF][START_REF] Almkvist | The method of differentiating under the integral sign[END_REF]; they reduce the resolution of the telescopic equation ( 25) to the computation of the rational solutions of a system of linear differential equations. The roots of this method can be traced back to Picard [START_REF] Picard | Sur les intégrales doubles de fonctions rationnelles dont tous les résidus sont nuls[END_REF] for n = 2. Algorithms from the third generation (3G) only use linear algebra, and are based on an idea that was first formulated by Apagodu and Zeilberger in [START_REF] Mohammed | Sharp upper bounds for the orders of the recurrences output by the Zeilberger and q-Zeilberger algorithms[END_REF][START_REF] Apagodu | Multi-variable Zeilberger and Almkvist-Zeilberger algorithms and the sharpening of Wilf-Zeilberger theory[END_REF], and has later been refined and generalized [START_REF] Koutschan | A fast approach to creative telescoping[END_REF][START_REF] Chen | Trading order for degree in creative telescoping[END_REF][START_REF] Chen | Order-degree curves for hypergeometric creative telescoping[END_REF][START_REF] Chen | A generalized Apagodu-Zeilberger algorithm[END_REF]. This approach is interesting not only because it is easier to implement and tends to run faster than earlier algorithms, but also because it is easy to analyze. A common drawback of these three generations of algorithms is that they all compute certificates, whose size is much bigger than that of telescopers. Moreover, 1G algorithms are slow, 2G algorithms have a bad or unknown complexity, and 3G algorithms do not necessarily output telescopers of minimal orders. However, already algorithms from the second generation are able to solve non-trivial problems.

Example 42 (Diagonal 3D Rook paths, cont.). Using notation from Example 34 and from the proof of Theorem 35, the aim is to construct a linear differential operator P(t, ∂ t ), and two rational functions R and S in Q(t, x, y) such that

P(G) = ∂R ∂x + ∂S ∂y .
Maple's implementation of Chyzak's algorithm is able to do this in a few seconds:

> G:=subs(y=y/x,z=t/y,1/(1-x/(1-x)-y/(1-y)-z/(1-z)))/y/x: > P,R,S:=op(op(Mgfun:-creative_telescoping(G,t::diff, [x::diff,y::diff])));

It outputs the differential equation P(∆) = 0 satisfied by ∆ = ∑ n≥0 D n t n , where

P = t(t -1)(64t -1)(3t -2)(6t + 1)∂ 3 t +(4608t 4 -6372t 3 + 813t 2 + 514t -4)∂ 2 t +4(576t 3 -801t 2 -108t + 74)∂ t ,
which helps proving a recurrence conjectured in [START_REF] Erickson | Enumerating rook and queen paths[END_REF].

4G Creative Telescoping. Algorithms from the fourth and most recent generation of creative telescoping algorithms are called reduction-based algorithms. Its roots are in works by Hermite [START_REF] Hermite | Sur l'intégration des fractions rationnelles[END_REF] and Picard [START_REF] Picard | Sur les intégrales doubles de fonctions rationnelles dont tous les résidus sont nuls[END_REF][START_REF] Picard | Théorie des fonctions algébriques de deux variables indépendantes[END_REF]. This approach was first applied to the integration of bivariate rational functions by Bostan, Chen, Chyzak and Li [START_REF] Bostan | Complexity of creative telescoping for bivariate rational functions[END_REF]. This first article generated a very active area of research [START_REF] Chen | Telescopers for rational and algebraic functions via residues[END_REF][START_REF] Bostan | Hermite reduction and creative telescoping for hyperexponential functions[END_REF][START_REF] Bostan | Creative telescoping for rational functions using the Griffiths-Dwork method[END_REF][START_REF] Chen | A modified Abramov-Petkovšek reduction and creative telescoping for hypergeometric terms[END_REF][START_REF] Huang | New bounds for hypergeometric creative telescoping[END_REF][START_REF] Lairez | Computing periods of rational integrals[END_REF][START_REF] Bostan | Efficient algorithms for mixed creative telescoping[END_REF][START_REF] Chen | Reduction-based creative telescoping for algebraic functions[END_REF][START_REF] Chen | Reduction-based creative telescoping for fuchsian D-finite functions[END_REF][START_REF] Bostan | Multiple binomial sums[END_REF].

Let us explain the principle of the method in the univariate case, that is when n = 1 in the telescopic Equation [START_REF] Banderier | Limit laws for basic parameters of lattice paths with unbounded jumps[END_REF].

The problem at hand is: given H = P/Q ∈ K(t, x), compute γ H(t, x)dx. The principle of the method originates from the Hermite reduction [START_REF] Hermite | Sur l'intégration des fractions rationnelles[END_REF], a procedure for computing a normal form of a univariate function modulo derivatives. Hermite introduced his method as a way to compute the algebraic part of the primitive of a univariate rational function without computing the roots of its denominator, as opposed to the classical partial fraction decomposition method.

By Hermite reduction, the integrand H can be written in reduced form

H = ∂ x (g) + a Q , where Q is the squarefree part of Q and deg x (a) < d := deg x (Q ).
The principle of the algorithm in [START_REF] Bostan | Complexity of creative telescoping for bivariate rational functions[END_REF] is then the following: 1. For i = 0, 1, . . . , d compute the Hermite reduction of ∂ i t (H):

∂ i t (H) = ∂ x (g i ) + a i Q , deg x (a i ) < d .
2. Find the first linear relation over Q(t) of the form ∑ r k=0 c k a k = 0. Then L = ∑ r k=0 c k ∂ k t is a telescoper, and ∑ r k=0 c k g k the corresponding certificate. The method has been extended to the multivariate case of periods of rationals integrals by Bostan, Lairez and Salvy [START_REF] Bostan | Creative telescoping for rational functions using the Griffiths-Dwork method[END_REF]. They have obtained the following result. Theorem 43 ([87]). Let H = P Q be a rational function in t and x = x 1 , . . . , x n and denote by d x the degree of Q w.r.t. x, and d t = max(deg t P, deg t Q). Assume

deg x P + n + 1 d x .
Then a telescoper for H can be computed using O(e 3n d 8n x d t ) operations in Q, uniformly in all the parameters. The minimal telescoper has order ≤ d n

x and degree O(e n d 3n x d t ). These size bounds are generically reached. There are three main ideas behind the proof of Theorem 43:

• in the generic case, a multivariate generalization of Hermite's reduction is used; it called the Griffiths-Dwork method [165, §3], [166, §8], [START_REF] Griffiths | On the periods of certain rational integrals. I, II[END_REF]; • in the general case, a deformation technique is used to reduce to the generic case, by an input perturbation using a new free variable; • fast linear algebra algorithms for polynomial matrices [START_REF] Storjohann | High-order lifting and integrality certification[END_REF]396] is used to deal with Macaulay matrices that encode Gröbner bases computations. The algorithm behind Theorem 43 is the first algorithm for creative telescoping with polynomial complexity in the generic size of the output Picard-Fuchs equation. It avoids the costly computation of certificates. This is crucial since, generically, certificates have size Ω(d n 2 /2 x ). Previous algorithms have (at least) doubly-exponential complexity, inherited from the fact that they need to compute certificates. A recent, and highly non-trivial, extension of the results in [START_REF] Bostan | Creative telescoping for rational functions using the Griffiths-Dwork method[END_REF] was given by Lairez [START_REF] Lairez | Computing periods of rational integrals[END_REF]. It tremendously improves the practical efficiency of the algorithm in [START_REF] Bostan | Creative telescoping for rational functions using the Griffiths-Dwork method[END_REF].

∑ k=0 4 k ( 2k k ) , n ∑ k=0 k ∑ j=0 n j 3 or n ∑ i=0 n ∑ j=0 i + j j 2 4n -2i -2j 2n -2i .
Many multiple sums can be cast into problems of rational integration by passing to generating functions. This observation was intensively used by Egorychev in his book [START_REF] Egorychev | Integral representation and the computation of combinatorial sums[END_REF], but its algorithmic consequences were studied only quite recently by Bostan, Lairez and Salvy [START_REF] Bostan | Multiple binomial sums[END_REF]. They defined a class of multi-indexed sequences called (multiple) binomial sums, which is closed under partial summation, and contains most of the sequences obtained by multiple summation of products of binomial coefficients and also all the sequences with algebraic generating function. Not every sum that creative telescoping can handle is a binomial sum: for example, among the three sums in Eq. ( 26), the second one and the third one are binomial sums but the first one is not, since it contains the inverse of a binomial coefficient. Yet many sums coming from combinatorics and number theory are binomial sums. The starting point is that integral representations of the generating function of a binomial sum can be computed in an automated way. The outcome is twofold. Firstly, the generating functions of univariate binomial sums are exactly the diagonals of rational power series; this equivalence characterizes binomial sums in an intrinsic way. All the theory of diagonals transfers to univariate binomial sums and gives many interesting arithmetic properties. Secondly, integral representations can be used to actually compute with binomial sums (e.g. find recurrence relations or prove identities automatically) via the computation of Picard-Fuchs equations.

Example 44. (A particular instance of Dixon's identity) We will simply illustrate the main points of the method in [START_REF] Bostan | Multiple binomial sums[END_REF] on the identity [START_REF] Banderier | The kernel method for lattice paths below a line of rational slope[END_REF] 2n

∑ k=0 (-1) k 2n k 3 = (-1) n (3n)! n! 3 .
The strategy is as follows: find an integral representation of the generating function of the left-hand side; simplify this integral representation using partial integration; use the simplified integral representation to compute a differential equation of which the generating function is solution; transform this equation into a recurrence relation; solve this recurrence relation. First of all, the binomial coefficient

( n k ) is the coefficient of x k in (1 + x) n . Cauchy's integral formula ensures that n k = 1 2πi γ (1 + x) n x k dx x ,
where γ is the circle x ∈ C |x| = 1 2 . Therefore, the cube of a binomial coefficient can be represented as a triple integral

2n k 3 = 1 (2πi) 3 γ×γ×γ (1 + x 1 ) 2n x k 1 (1 + x 2 ) 2n x k 2 (1 + x 3 ) 2n x k 3 dx 1 x 1 dx 2 x 2 dx 3 x 3 .
As a result, the generating function y(t) of the left-hand side of Equation ( 27) equals

= 1 (2iπ) 3 γ 3 ∞ ∑ n=0 2n ∑ k=0 t 3 ∏ i=1 (1 + x i ) 2 n -1 x 1 x 2 x 3 k dx 1 x 1 dx 2 x 2 dx 3 x 3 = 1 (2iπ) 3 γ 3 ∞ ∑ n=0 t 3 ∏ i=1 (1 + x i ) 2 n 1 - -1 x 1 x 2 x 3 2n+1 1 + 1 x 1 x 2 x 3 dx 1 x 1 dx 2 x 2 dx 3 x 3 = 1 (2iπ) 3 γ 3 x 1 x 2 x 3 -t ∏ 3 i=1 (1 + x i ) 2 dx 1 dx 2 dx 3 x 2 1 x 2 2 x 2 3 -t ∏ 3 i=1 (1 + x i ) 2 1 -t ∏ 3 i=1 (1 + x i ) 2
.

The partial integral with respect to x 3 along the circle |x 3 | = 1 2 is the sum of the residues of the rational function being integrated at the poles whose modulus is less than 1 2 . When |t| is small and |x 1 | = |x 2 | = 1 2 , the poles coming from the factor x 2 1 x 2 2 x 2 3t ∏ 3 i=1 (1 + x i ) 2 all have a modulus that is smaller than 1 2 : they are asymptotically proportional to |t| 1/2 . In contrast, the poles coming from the factor 1t ∏ 3 i=1 (1 + x i ) 2 behave like |t| -1/2 and have all a modulus that is bigger than 1 2 . In particular, any two poles that come from the same factor are either both asymptotically small or both asymptotically large. This implies that the partial integral is a rational function of t, x 1 and x 2 ; and we compute that

y(t) = 1 (2iπ) 2 γ×γ x 1 x 2 dx 1 dx 2 x 2 1 x 2 2 -t(1 + x 1 ) 2 (1 + x 2 ) 2 (1 -x 1 x 2 )
2 . This formula echoes the original proof of [START_REF] Dixon | On the sum of the cubes of the coefficients in a certain expansion by the binomial theorem[END_REF] in which the left-hand side of ( 27) is expressed as the coefficient of (xy) 4n in ((1y 2 )(1z 2 )(1y 2 z 2 )) 2n . Using any algorithm described in §3.2.3 that performs definite integration of rational functions reveals a differential equation satisfied by y(t): t(27t + 1)y + (54t + 1)y + 6y = 0.

Looking at the coefficient of t n in this equality leads to the recurrence relation 3(3n + 2)(3n + 1)u n + (n + 1) 2 u n+1 = 0, where u n = ∑ 2n k=0 (-1) k ( 2n k ) 3 . Since u 0 = 1, it proofs Dixon's identity [START_REF] Banderier | The kernel method for lattice paths below a line of rational slope[END_REF].

Note that the method avoids the computation of certificates; this nice feature is inherited from the computation of Picard-Fuchs equations for periods of rational integrals, which can be achieved efficiently without computing the corresponding certificate and without introducing spurious singularities ( §3.2.3). This should be contrasted with the usual creative telescoping methods for sums ( §3.2.1).

Creative Telescoping for quarter plane walks.

Let us now turn back to quarter plane walks with small steps. We focus on models 1-19 in Fig. 7, and to Theorem 14. We write F(t; x, y) for the full generating function F S (t; x, y), where S is one of the 19 models.

Using the kernel method, Bousquet-Mélou and Mishna showed in [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF]Prop. 8] that the generating function F(t; x, y) can be written in the form where the Hadamard product denoted x,y is the term-wise product of two series, while the diagonal operator Diag x,x selects those terms with equal exponents of x and x . This argument also implies an algorithm for computing linear differential equations satisfied by F(t; x, y), since diagonals can be computed using creative telescoping. Therefore, from (28) one could, in principle, determine differential equations for F(t; x, y). However, the direct use of ( 29) in our context leads to infeasible computations; worse, the intermediate algebraic objects involved in the calculations would probably have too large sizes to be merely written and stored. This is really unfortunate, since our need is mere evaluations of the diagonals in [START_REF] Beckermann | A uniform approach for the fast computation of matrix-type Padé approximants[END_REF] at specific values for x and y.

Example 45. (King Walks in the Quarter Plane) We illustrate the approach on the king walks (model 4 with S = in Fig. 7). The first terms of the length generating function F(t; 1, 1) read (see http://oeis.org/A151331) F(t; 1, 1) = 1 + 3t + 18t 2 + 105t 3 + 684t 4 + 4550t 5 + 31340t 6 + 219555t 7 + • • • , and we describe the method used in [START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF] to obtain the closed formula (9) for it.

First, the kernel equation ( 7) reads [START_REF] Benoit | Quasi-optimal multiplication of linear differential operators[END_REF] xyJ(x, y)F(x, y) = xytx(x + 1 + x)F(x, 0)ty(y + 1 + ȳ)F(0, y) + tF(0, 0), where F(x, y) ≡ F(t; x, y), x := 1/x, ȳ := 1/y and J(x, y) is the Laurent polynomial J(x, y) = 1t ∑ (i,j)∈S

x i y j = 1t(xy + y + xy + x + x + x ȳ + ȳ + x ȳ).

The group of S has order 4: it contains the elements (x, y), ( x, y), ( x, ȳ), (x, ȳ), which leave invariant J(t; x, y). Applying these rational transformations to the kernel equation (30) yields the four relations:

xyJ(x, y)F(x, y) = xytx(x + 1 + x)F(x, 0)ty(y + 1 + ȳ)F(0, y) + tF(0, 0), -xyJ(x, y)F( x, y) =xy + t x(x + 1 + x)F( x, 0) + ty(y + 1 + ȳ)F(0, y) -tF(0, 0),

x ȳJ(x, y)F( x, ȳ) = x ȳt x(x + 1 + x)F( x, 0)t ȳ(y + 1 + ȳ)F(0, ȳ) + tF(0, 0), -xyJ(x, y)F( x, y) = -x ȳ + tx(x + 1 + x)F(x, 0) + t ȳ(y + 1 + ȳ)F(0, ȳ) -tF(0, 0). Upon adding up these equations, all terms in the right-hand side involving F disappear, resulting in xyF(x, y) -xyF( x, y) + x ȳF( x, ȳ)x ȳF(x, ȳ) = J(x, y) -1 (xyxy + x ȳx ȳ) . Now, the main observation is that on the left-hand side, all terms except the first one involve negative powers either of x or of y. Therefore, extracting positive parts expresses the generating series xyF(x, y) as the positive part (w.r.t. x and y) of a trivariate rational function: [START_REF] Bernardi | Bijective counting of Kreweras walks and loopless triangulations[END_REF] xyF(x, y) = [x > ][y > ] xyxy + x ȳx ȳ 1t(xy + y + y x + x + x ȳ + ȳ + x ȳ + x) .

Up to this point, the reasoning is borrowed from Bousquet-Mélou's and Mishna's article [START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF]. Combined with Theorem 33, it already implies that F(x, y) is D-finite; in particular, F(1, 1) is also D-finite. Our aim is to refine this qualitative result, and explicitly obtain a linear differential equation satisfied by F(1, 1). Starting from [START_REF] Bernardi | Bijective counting of Kreweras walks and loopless triangulations[END_REF] and following more closely Lipshitz' encoding [START_REF] Lipshitz | The diagonal of a D-finite power series is D-finite[END_REF], a first observation is that F(x, y) is equal to the iterated diagonal Diag x 1 ,x 2 Diag y 1 ,y 2 of the rational function [START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF] x 2 y 2 (x 1 y 1 -x1 y 1 + x1 ȳ1x 1 ȳ1 ) (1x 2 )(1y 2 )(1t(x 1 y 1 + y 1 + y 1 x1 + x1 + x1 ȳ1 + ȳ1 + x 1 ȳ1 + x 1 )) .

However, this computation is too difficult, and exceeds by far the limits of the best existing algorithms for diagonals. The reason is that differential equations w.r.t. t and with polynomial coefficients in x, y, t are really huge, so the main limitation of algorithms computing (32) already comes from the size of the output. Another weakness of the diagonal encoding [START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF] is that it does not provide direct access to the univariate series F(1, 1), since taking diagonals and specializing variables are operations that do not commute.

To circumvent these difficulties and to make the computation feasible, the key idea in [START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF] is to encode the positive part in [START_REF] Bernardi | Bijective counting of Kreweras walks and loopless triangulations[END_REF] as a formal residue: [START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF] F(α, β) = [x -1 y -1 ] xyxy + x ȳx ȳ (1αx)(1βy)(1t(xy + y + y x + x + x ȳ + ȳ + x ȳ + x))

.

The formal proof of this encoding is delicate. The advantage of (33) over ( 32) is twofold. On the one hand, the residue computation can be carried out by using a single call to the creative-telescoping algorithm for rational functions, while the diagonal computation [START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF] has two steps, the first for a rational function in five variables, the second for an algebraic function in four variables. On the other hand, and more importantly, taking residues commutes with specialization, contrarily to positive parts and diagonals. Therefore, the generating series for walks F(1, 1) is

F(1, 1) = [x -1 y -1 ] xy -xy + x ȳ -x ȳ (1 -x)(1 -y)(1 -t(xy + y + y x + x + x ȳ + ȳ + x ȳ + x))
, and a differential equation L(F(1, 1)) = 0 can now be computed by creative telescoping:

(34) L = t 2 (1 + 4t)(8t -1)(2t -1)(1 + t)∂ [START_REF] Kunwar | Second order differential equations with hypergeometric solutions of degree three[END_REF][START_REF] Imamoglu | Computing hypergeometric solutions of second order linear differential equations using quotients of formal solutions[END_REF][START_REF] Imamoglu | Computing hypergeometric solutions of second order linear differential equations using quotients of formal solutions and integral bases[END_REF]) allow to express f (t) as

f (t) = 1 t(1 + 4t) 3 • 2 F 1 3 2 3 2 2 16t(1 + t) (1 + 4t) 2 .
Finally, solving the equation d/dt F(1, 1) + F(1, 1)/t = f (t) yields formula [START_REF] Almkvist | The method of differentiating under the integral sign[END_REF]. Similarly, for indeterminates α and β we obtain the formal residue representations for F(α, 0) and F(0, β), and creative-telescoping techniques still allow the effective computation of differential operators for F(α, 0), resp. for F(0, β). Owing to the additional symbolic indeterminate, the computations are much harder than for F(1, 1), but still feasible. Each of the resulting differential operators factors again, this time as a product of an order-two operator and of three order-one operators. Moreover, the second-order operators are again solvable in terms of 2 F 1 functions. Finally, a closed formula for F(α, β) is obtained from the closed formulas for F(α, 0) and F(0, β) via the kernel equation [START_REF] Benoit | Quasi-optimal multiplication of linear differential operators[END_REF]. This detour is computationally crucial, since performing creative telescoping directly on the five-variable rational function from [START_REF] Bernardi | Counting quadrant walks via Tutte's invariant method[END_REF] is not feasible even using today's best algorithms.

A similar reasoning applies to any of the 19 models in Fig. 7 with finite group and non-zero orbit sum, and this allows to prove Theorem 14 with the help of the fundamental equation GF = PositivePart orbit sum kernel .

3.5. Back to the exercise in §1.1. To conclude, we come back to the problem stated at the very beginning of the memoir, for which we have guessed the answer in §2.6. Recall that S denotes the step set {↑, ←, }. For convenience, we will continue to use the shortcut notation x = 1/x, ȳ = 1/y. It satisfies the kernel equation [START_REF] Albert | Permutations sortable by two stacks in parallel and quarter plane walks[END_REF], which reads: [START_REF] Bertrand | Solution d'un problème[END_REF] 1t(y + x + x ȳ) xyQ(x, y) = xytx 2 Q(x, 0) -tyQ(0, y).

We are interested in the generating function of diagonal returns B(t) = [x 0 ] Q(x, x). It satisfies a functional equation very similar to [START_REF] Bertrand | Solution d'un problème[END_REF], namely [START_REF] Beukers | Padé-approximations in number theory[END_REF] 1t(y + x + x ȳ) xyH(x, y) = xytx 2 H(x, 0).

A functional equation for

This time, we are interested in A(t) = [x 0 ] H(x, 0), the generating function of excursions in the upper half-plane.

3.5.3. The kernel method for Z × N. We solve Eq. 37 by using the same technique as we did for Dyck walks (Equation (3) from Example 5).

Let

y 0 = x -t -(t -x) 2 -4t 2 x 3 2tx = xt + t 2 + (x 2 + x)t 3 + (3x + x2 )t 4 + • • • be the (unique) root in Q[x, x][[t]
] of K(x, y 0 ) = 0, where K(x, y) = 1t(y + x + x ȳ). Then plugging y 0 in (37) yields 0 = K(x, y 0 )yH(x, y 0 ) = y 0 -txH(x, 0), and thus H(x, 0) = y 0 tx and A(t) = x 0 y 0 tx . This allows to express A(t) as a period of an algebraic integral. A differential equation satisfied by A(t) can then be computed using creative telescoping: > y0:= -sqrt((t-x)^2 -4*t^2*x^3)/(2*t*x): > DEtools [Zeilberger](1/x * y0/(t*x), t, x, Dt) [START_REF] Abel | OEuvres complètes. Tome II[END_REF]; which proves the equation (27t 4t)A (t) + (108t 3 -4)A (t) + 54t 2 A(t) = 0, or equivalently, the recurrence relation on its coefficients: 27(n + 2)(n + 1)a n = (n + 6)(n + 3)a n+3 . Then, applying the method in [START_REF] Bostan | Hypergeometric expressions for generating functions of walks with small steps in the quarter plane[END_REF] allows to express B(t) as a residue:

B(t) = [x 0 ]Q(x, x) = [u -1 v -1 z -1 ] ū v -u v2 + u 2 v -uv + ūv 2 -ū2 v z(1 -zu)(1 -v z)(1 -t( v + u + ūv))
.

Finally, multivariate Creative Telescoping proves a differential equation for B(t):

> OS := x*y -y^2/x + y/x^2 -1/x/y + x/y^2 -x^2/y; > ker := 1-t*(y + 1/x + x/y); > S:=normal(subs({x=1/u,y=1/v}, OS/ker)/(1-z*u)/(1-v/z)/z); > Mgfun:-creative_telescoping(S,t::diff,[z::diff,u::diff,v::diff]):

namely (27t 4t)B (t) + (108t 3 -4)B (t) + 54t 2 B(t) = 0.

Conclusion.

We have proved that A(t) and B(t) are both solutions of (27t 4t)y (t) + (108t 3 -4)y (t) + 54t 2 y(t) = 0.

Solving this equation in closed form proves that

A(t) = B(t) = 2 F 1 1/3 2/3 2 27 t 3 = ∞ ∑ n=0 (3n)! n! 3 t 3n n + 1 .
Thus the two sequences are equal to

a 3n = b 3n = (3n)! n! 2 • (n + 1)!
, and a m = b m = 0 if 3 does not divide m.

Figure 1 .

 1 Figure 1. Some discarded models: trivial; intrinsic to the half-plane; symmetric.

Figure 2 .

 2 Figure 2. The 79 models of small-step walks in the quarter plane: 74 non-sigular, 5 singular.

Figure 3 .

 3 Figure 3. The most basic classes of power series, and their dependencies.

Figure 4 .

 4 Figure 4. The Beukers-Heckman interlacing criterion [40] at work.

∑

Figure 5 .

 5 Figure 5. The simple walk in the cones with angle 45 • and 135 • : Gouyou-Beauchamps and Gessel walks.

1. 14 .

 14 Algebraic reformulation: solving a functional equation. Gessel's problem admits the following purely algebraic reformulation, which should be seen as a quarter-plane analogue of Equation (3) from Example 5. If G(t; x, y) ∈ Q[x, y][[t]] denotes the full generating function for Gessel walks in the quarter plane then a simple inclusion-exclusion reasoning represented pictorially in Fig. 6 implies that G(t; x, y) satisfies a functional equation called the kernel equation G (t;

Figure 6 .

 6 Figure 6. The functional equation for Gessel walks in the quarter plane, pictorially.

Figure 9 .

 9 Figure 9. Examples of models with groups of orders 4, 6, 8 and ∞, respectively.

Fig. 7

 7 admit full generating functions that are algebraic. Moreover, closed formulas exist for them. For the three models 20-22 related to the Kreweras model, such formulas are displayed in[101, §6]. The most difficult case among these four is model 23(Gessel's), for which Theorem 13 provides a closed-form expression.

Figure 11 .

 11 Figure 11. Rotations of a scarecrow: models with zero drift that have a non-D-finite generating function.

Figure 12 .

 12 Figure 12. The 9 models with a non-D-finite but D-algebraic generating function.

|G S |<∞: 23 nonzero 23 .

 2323 Extensions and open questions.

  |G S | = ∞: 20 634[START_REF] Du | Infinite orders and non-D-finite property of 3-dimensional lattice walks[END_REF] Th. 1.3]

  Figure 14. Guessing differential operators for G(t; a, 0), for prime p and a ∈ F p : minimal-order operator (blue point above the hyperbola) obtained as gcrd of several non-minimal operators (blue points below the hyperbola). Points below the hyperbola correspond to operators obtainable by Hermite-Padé approximation with 1000 terms.

2 .

 2 (Uniqueness) By the implicit function theorem, this polynomial P admits a root r(t) ∈ Q[[t]] with r(0) = 1. Since P(T, 0) = T -1 has a single root in C, the series r(t) is the unique root of P in C[[t]].

Figure 15 .

 15 Figure 15. Algebraicity of Kreweras walks: a computerized proof in a nutshell.

2 . 5 .

 25 . (D-finiteness) Discover and certify a differential equation L for f (t) of order 11 and degree 73 2. (Local analysis) L is Fuchsian and has a logarithmic singularity at t = 0 3. (Bounds) If ord(L f min ) ≤ 10, then L f min has coefficients of degrees ≤ 580 4. (Guessing) Differential Hermite-Padé approximants rule out this possibility 5. (Conclusion) Thus, L f min = L, and so f is transcendental. The bounds in Step 3 are the mathematical heart of the proof: they are obtained by using the Fuchsianity of L, and by bounding the apparent singularities of factors of L via Fuchs' equality, cf. [379, §4.4.1] and [334, §20]. Inside the toolbox: Hermite-Padé approximants. We now have a quick closer look at what is hidden behind guessing: Hermite-Padé approximants.

Figure 17 .

 17 Figure 17. The diagonal of a bivariate power series (on the left) viewed as a residue (on the right).

1, 2 ,

 2 14, 106, 838, 6802, 56190, 470010, . . .

Figure 18 .

 18 Figure 18. Number d N of diagonal 2D Rook paths from lower-left corner to N × N upper-right corner.

3. 2 . 1 .

 21 Creative Telescoping for sums. Let us explain the basic principle of the method on the simplest example possible. Denote by I n the definite sum I n :=

>

  pN:=1/4: pS:=1/4: pE:=1/4 + e: pW:=1/4 -e: > U:=binomial(2*n,2*k)* binomial(2*k,k)*pE^k*pW^k * binomial(2*n-2*k,n-k)*pN^(n-k)*pS^(n-k): > SumTools[Hypergeometric][Zeilberger](U,n,k,Sn): > collect(%[1], Sn, factor);

( 28 )

 28 F(t; x, y) = 1 xy [x > ][y > ] N(x, y) 1 -tS(x, y) COMPUTER ALGEBRA FOR LATTICE PATH COMBINATORICS 59where N(x, y) and S(x, y) are certain Laurent polynomials in y with coefficients that are rational functions in x. The intended reading of[START_REF] Barton | Some aspects of the random sequence[END_REF] is: first interpret N(x, y)/(1 -tS(x, y)) as an element of Q(x)[y, 1/y][[t]]; let [y > ] act term by term, obtaining a series in Q(x)[y][[t]] that actually belongs to Q[x, 1/x][y][[t]] for all cases in Figure7; then let [x > ] act term by term, finally obtaining an element ofQ[x][y][[t]].In this reading, the composition [x > ][y > ] of positive-part operators is only applied to Laurent polynomials, for which it is well-defined, in a unique way.As pointed out by Bousquet-Mélou and Mishna, Equation (28) already implies the D-finiteness of F(t; x, y), by Theorem 33 and since positive parts can be encoded as diagonals. To be more specific, the positive part[x > ][y > ]R(t; x, y) of a formal power series R ∈ Q[[x,y, t]] can be encoded as (29) x 1x y 1y x,y R(t; x, y) = Diag x,x Diag y,y x 1x y 1y R(t; x , y ),

3. 5 . 1 . 2 Q

 512 A functional equation for S-walks in N 2 . Let us consider the full generating function for S-walks in N i,j t n x i y j ∈ Q[x, y][[t]].

  S-walks in Z × N. Similarly, let H(t; x, y) ≡ H(x, y) denote the full generating function for S-walks in Z × N, ; i, j)t n x i y j ∈ Q[x, x, y][[t]].

3. 5 . 4 .

 54 The kernel method for N 2 . The inventory χ(x, y) = x + y + x ȳ of S is left unchanged by the involutions Φ : (x, y) → ( xy, y) and Ψ : (x, y) → (x, x ȳ) .

  which generate a finite dihedral group D 3 of order 6:( xy, y) (x, x ȳ) ( xy, x) ( ȳ, x ȳ) x)Letting the group act on the kernel equation[START_REF] Bertrand | Solution d'un problème[END_REF] gives six equations, whose alternate sum gives birth to the orbit equation:xyQ(x, y)xy 2 Q( xy, y) + x2 yQ( xy, x) x ȳQ( ȳ, x) + x ȳ2 Q( ȳ, x ȳ)x 2 ȳQ(x, x ȳ) = xyxy 2 + x2 yx ȳ + x ȳ2x 2 ȳ 1t(y + x + x ȳ)Extracting the part with positive powers of x and y like in(3.4) givesxyQ(x, y) = [x >0 y >0 ] xyxy 2 + x2 yx ȳ + x ȳ2x 2 ȳ 1t(y + x + x ȳ) .

  and c, either {ra} ≤ {rc} < {rb} or {rb} ≤ {rc} < {ra}, where {x} denotes the fractional part xx of x. Let {a 1 , . . . , a k } and {b 1

			obtained an alternative
	arithmetic proof of Schwarz' result, which is more elementary and algorithmic.
	Assume w.l.o.g. that a, b, c ∈ Q such that a, b, c -a, c -b / ∈ Z. Then Errera's
	criterion states that 2 F 1	a b c	t is algebraic if and only if for every r coprime with
	the denominators of a, b For instance, this allows to prove
	immediately that
	• • 2 F 1 2 F 1 and that -1 2 -1 6 2 3 2t 1 12 5 12 1 1728 t = 1 + 60 t + 39780 t 2 + 38454000 t 3 + • • • is transcen-16 t -1 = 1 + 2 t + 11 t 2 + 85 t 3 + 782 t 4 + • • • is algebraic,
	dental.
	A generalization of this result, which completely solves Fuchs' question, was ob-
	tained by Beukers and Heckman in 1989 [40].
	Theorem 8.

  x, y) =1 + t xy + x +

	1 xy	+	1 x	G(t; x, y)

  , F(t; x, y) was proved transcendental byMishna [308, Th. 2.5]; for model 17, Mishna[308, §2.3.3] and Bousquet-Mélou and Mishna [101, §5.2], showed that F(t; x, y) and F(t; 0, 0) are transcendental and that F(t; 1, 1) is algebraic; for model 18, F(t; 1, 1) was proved algebraic by Bousquet-Mélou and Mishna[101,

§5.2]; for model 19, Bousquet-Mélou and Mishna

[101, §5.3]

showed that F(t; 0, 0), F(t; 0, 1), F(t; 1, 0) and F(t; 1, 1) are transcendental. The first unified transcendence proof for F(t; x, y) applying to all 19 models is by Fayolle and Raschel [181, Theorem 1.1], although they attribute that result to Bousquet-Mélou and Mishna

[START_REF] Bousquet-Mélou | Walks with small steps in the quarter plane[END_REF]

.

  and for the Weyl algebra F p [t] ∂ t [217, 376, 55, 71, 30, 75, 377]. They are not needed for simple examples, but they become of crucial help in the treatment of examples of critical sizes, such as for the computations involved in Gessel's model, see Example 24.

  ).

	> gfun:-listtoalgeq([1,1,2,5,14], y(t));
	2
	[1 -y(t) + t y(t) , ogf]
	> gfun:-listtodiffeq([1,1,2,5,14,42,132,429], y(t), [ogf]);
	FAIL
	> gfun:-listtodiffeq([1,1,2,5,14,42,132,429,1430,4862], y(t), [ogf

  3 t 3 + 30 t 6 + 420 t 9 + 6930 t 12 + 126126 t 15 + 2450448 t 18 +49884120 t 21 + 1051723530 t 24 + 22787343150 t 27 + O(t 30 ).

	(19)

  30 t 6 + 420 t 9 + 6930 t 12 + 126126 t 15 + 2450448 t 18 +49884120 t 21 + 1051723530 t 24 + 22787343150 t 27 + O(t 30 ).

  As explained in §3.2.1, creative telescoping allows to prove identities like Dixon's (first item in Example 36), and to deal with definite sums like

	COMPUTER ALGEBRA FOR LATTICE PATH COMBINATORICS	57
	3.3. Binomial sums. (26)

n

  3 t + t(200t 3 + 576t 4 -33t -252t 2 + 5)∂ 2 t + 4(22t 3 -117t 2 -12t + 288t 4 + 1)∂ t + 384t 3 -12 -144t -72t 2 .COMPUTER ALGEBRA FOR LATTICE PATH COMBINATORICS 61 Note that this is precisely the differential operator guessed in[START_REF] Bostan | Automatic classification of restricted lattice walks[END_REF].Moreover, factorization algorithms for linear differential operators[START_REF] Grigor Ev | Complexity of factoring and calculating the GCD of linear ordinary differential operators[END_REF][START_REF] Singer | Testing reducibility of linear differential operators: a group-theoretic perspective[END_REF][START_REF] Bronstein | An introduction to pseudo-linear algebra[END_REF][START_REF] Van Hoeij | Factorization of differential operators with rational functions coefficients[END_REF] can be used to prove that L = L 2 L 1 , whereL 1 = ∂ t + 1/t and (35) L 2 = t 2 (1 + 4t)(1 -8t)(1 -2t)(1 + t)∂ 2 t + 2t(256t 4 + 80t 3 -111t 2 -14t + 2)∂ t + 768t 4 + 8t 3 -306t 2 -30t + 2. -1 + 6 + 54t + 420t 2 + 3420t 3 + 27300t 4 + O(t 6 )is a solution of L 2 . Starting from the second order operator L 2 , algorithmic methods explained in[77, §2.6] (see also

	It follows that the Laurent power series
	f (t) =	dF dt	(1, 1) +	F(1, 1) t	= t

§ Bousquet-Mélou and Mishna[101, §3] do so for the 51 non-singular models, but F. Chyzak [private communication] points out that the argument still works on some iterate of θ.

¶ In all that follows, we have used the version 3.76 (July 2015) of the package gfun, available at http://perso.ens-lyon.fr/bruno.salvy/software/the-gfun-package/.

* * * By highly optimizing this argument[START_REF] Bostan | Explicit formula for the generating series of diagonal 3D rook paths[END_REF] reduces the problem to a kernel computation of a polynomial matrix of size 8917 × 9139, with entries in Q[x] of degree at most 37: these sizes are still too high to be dealt with in practice.

† † † ISSAC, the International Symposium on Symbolic and Algebraic Computation, is the premier conference for research in symbolic computation and computer algebra http://www.issac-symposium.org.