Discrete Differential Equations (DDEs) are functional equations that relate algebraically a power series F(t, u) in t with polynomial coefficients in a “catalytic” variable u and the specializations, say at u = 1, of F(t, u) and of some of its partial derivatives in u. If a DDE is of a fixed-point type then its solution F(t, u) is unique, and an elegant result by Bousquet-Mélou and Jehanne implies that F(t, u) is an algebraic power series. Last year, Bostan et al. initiated a systematic algorithmic study of DDEs of order 1. We generalize this study to DDEs of arbitrary order. First, we propose nontrivial extensions of algorithms based on polynomial elimination and on the guess-and-prove paradigm. Second, we design two brand-new algorithms that exploit the special structure of the underlying polynomial systems. Last, but not least, we report on implementations that are able to solve highly challenging DDEs with a combinatorial origin.

The Nth power of a polynomial matrix of fixed size and degree can be computed by binary powering as fast as multiplying two polynomials of linear degree in N. When Fast Fourier Transform (FFT) is available, the resulting complexity is softly linear in N, i.e. linear in N with extra logarithmic factors. We show that it is possible to beat binary powering, by an algorithm whose complexity is purely linear in N, even in absence of FFT. The key result making this improvement possible is that the entries of the Nth power of a polynomial matrix satisfy linear differential equations with polynomial coefficients whose orders and degrees are independent of N. Similar algorithms are proposed for two related problems: computing the Nth term of a C-finite sequence of polynomials, and modular exponentiation to the power N for bivariate polynomials.

Discrete differential equations of order 1 relate polynomially a power series F(t,u) in t with polynomial coefficients in a ''catalytic'' variable~u and one of its specializations, say F(t,u). Such equations are ubiquitous in combinatorics, notably in the enumeration of maps and walks. When the solution F is unique, a celebrated result by Bousquet-Mélou and Jehanne, reminiscent of Popescu's theorem in commutative algebra, states that F is algebraic. We address algorithmic and complexity questions related to this result. In generic situations, we first revisit and analyze known algorithms, based either on polynomial elimination or on the guess-and-prove paradigm. We then design two new algorithms: the first has a geometric flavor, the second blends elimination and guess-and-prove. In the general case (no genericity assumptions), we prove that the total arithmetic size of the algebraic equations for $F(t,1)$ is bounded polynomially in the size of the input discrete differential equation, and that one can compute such equations in polynomial time.

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 (e.g., algebraicity versus transcendence) and on explicit formulas for generating functions of walks with small steps in the quarter plane. In doing so, we emphasize the algorithmic nature of the methodology, especially two important paradigms: "guess-and-prove" and "creative telescoping".

We investigate the connection between properties of formal languages and properties of their generating series, with a focus on the class of holonomic power series. We first prove a strong version of a conjecture by Castiglione and Massazza: weakly-unambiguous Parikh automata are equivalent to unambiguous two-way reversal bounded counter machines, and their multivariate generating series are holonomic. We then show that the converse is not true: we construct a language whose generating series is algebraic (thus holonomic), but which is inherently weakly-ambiguous as a Parikh automata language. Finally, we prove an effective decidability result for the inclusion problem for weakly-unambiguous Parikh automata, and provide an upper-bound on its complexity.

In 1977, Strassen invented a famous baby-step / giant-step algorithm that computes the factorial N! in arithmetic complexity quasi-linear in [EQUATION]. In 1988, the Chudnovsky brothers generalized Strassen's algorithm to the computation of the N-th term of any holonomic sequence in the same arithmetic complexity. We design q-analogues of these algorithms. We first extend Strassen's algorithm to the computation of the q-factorial of N, then Chudnovskys' algorithm to the computation of the N-th term of any q-holonomic sequence. Both algorithms work in arithmetic complexity quasi-linear in [EQUATION]. We describe various algorithmic consequences, including the acceleration of polynomial and rational solving of linear q-differential equations, and the fast evaluation of large classes of polynomials, including a family recently considered by Nogneng and Schost.

Hermite reduction is a classical algorithmic tool in symbolic integration. It is used to decompose a given rational function as a sum of a function with simple poles and the derivative of another rational function. We extend Hermite reduction to arbitrary linear differential operators instead of the pure derivative, and develop efficient algorithms for this reduction. We then apply the generalized Hermite reduction to the computation of linear operators satisfied by single definite integrals of D-finite functions of several continuous or discrete parameters. The resulting algorithm is a generalization of reduction-based methods for creative telescoping.

We address the question of computing one selected term of an algebraic power series. In characteristic zero, the best algorithm currently known for computing the~Nth coefficient of an algebraic series uses differential equations and has arithmetic complexity quasi-linear in √N. We show that over a prime field of positive characteristic p, the complexity can be lowered to O(log N). The mathematical basis for this dramatic improvement is a classical theorem stating that a formal power series with coefficients in a finite field is algebraic if and only if the sequence of its coefficients can be generated by an automaton. We revisit and enhance two constructive proofs of this result for finite prime fields. The first proof uses Mahler equations, whose sizes appear to be prohibitively large. The second proof relies on diagonals of rational functions; we turn it into an efficient algorithm, of complexity linear in log N and quasi-linear in p.

The p-curvature of a system of linear differential equations in positive characteristic p is a matrix that measures how far the system is from having a basis of polynomial solutions. We show that the similarity class of the p-curvature can be determined without computing the p-curvature itself. More precisely, we design an algorithm that computes the invariant factors of the p-curvature in time quasi-linear in √ p. This is much less than the size of the p-curvature, which is generally linear in p. The new algorithm allows to answer a question originating from the study of the Ising model in statistical physics.

Creative telescoping is a powerful computer algebra paradigm-initiated by Doron Zeilberger in the 90's- for dealing with definite integrals and sums with parameters. We address the mixed continuous-discrete case, and focus on the integration of bivariate hypergeometric-hyperexponential terms. We design a new creative telescoping algorithm operating on this class of inputs, based on a Hermite-like reduction procedure. The new algorithm has two nice features: it is efficient and it delivers, for a suitable representation of the input, a minimal-order telescoper. Its analysis reveals tight bounds on the sizes of the telescoper it produces.

We design an algorithm for computing the p-curvature of a differential system in positive characteristic p. For a system of dimension r with coefficients of degree at most d, its complexity is O~ (p d rω) operations in the ground field (where ω denotes the exponent of matrix multiplication), whereas the size of the output is about p d r2. Our algorithm is then quasi-optimal assuming that matrix multiplication is (i.e. ω = 2). The main theoretical input we are using is the existence of a well-suited ring of series with divided powers for which an analogue of the Cauchy--Lipschitz Theorem holds.

The diagonal of a multivariate power series $F$ is the univariate power series DiagF generated by the diagonal terms of F. Diagonals form an important class of power series; they occur frequently in number theory, theoretical physics and enumerative combinatorics. We study algorithmic questions related to diagonals in the case where F is the Taylor expansion of a bivariate rational function. It is classical that in this case DiagF is an algebraic function. We propose an algorithm that computes an annihilating polynomial for DiagF. Generically, it is its minimal polynomial and is obtained in time quasi-linear in its size. We show that this minimal polynomial has an exponential size with respect to the degree of the input rational function. We then address the related problem of enumerating directed lattice walks. The insight given by our study leads to a new method for expanding the generating power series of bridges, excursions and meanders. We show that their first N terms can be computed in quasi-linear complexity in N, without first computing a very large polynomial equation.

Let V ∈ Q(i)(a1,..., an)(q1, q2) be a rationally parametrized planar homogeneous potential of homogeneity degree k ≠ −2, 0, 2. We design an algorithm that computes polynomial necessary conditions on the parameters (a1,..., an) such that the dynamical system associated to the potential V is integrable. These conditions originate from those of the Morales-Ramis-Simó integrability criterion near all Darboux points. The implementation of the algorithm allows to treat applications that were out of reach before, for instance concerning the non-integrability of polynomial potentials up to degree 9. Another striking application is the first complete proof of the non-integrability of the collinear three body problem.

We discuss theoretical and algorithmic questions related to the p-curvature of differential operators in characteristic p. Given such an operator L, and denoting by Ξ(L) the characteristic polynomial of its p-curvature, we first prove a new, alternative, description of Ξ(L). This description turns out to be particularly well suited to the fast computation of Ξ(L) when p is large: based on it, we design a new algorithm for computing Ξ(L), whose cost with respect to p is Õ(p0.5) operations in the ground field. This is remarkable since, prior to this work, the fastest algorithms for this task, and even for the subtask of deciding nilpotency of the p-curvature, had merely slightly subquadratic complexity Õ(p1.79).

We present a new reduction algorithm that simultaneously extends Hermite's reduction for rational functions and the Hermite-like reduction for hyperexponential functions. It yields a unique additive decomposition that allows to decide hyperexponential integrability. Based on this reduction algorithm, we design a new algorithm to compute minimal telescopers for bivariate hyperexponential functions. One of its main features is that it can avoid the costly computation of certificates. Its implementation outperforms Maple's function DEtools[Zeilberger]. We also derive an order bound on minimal telescopers that is tighter than the known ones.

Uncoupling algorithms transform a linear differential system of first order into one or several scalar differential equations. We examine two approaches to uncoupling: the cyclic-vector method (CVM) and the Danilevski-Barkatou-Zürcher algorithm (DBZ). We give tight size bounds on the scalar equations produced by CVM, and design a fast variant of CVM whose complexity is quasi-optimal with respect to the output size. We exhibit a strong structural link between CVM and DBZ enabling to show that, in the generic case, DBZ has polynomial complexity and that it produces a single equation, strongly related to the output of CVM. We prove that algorithm CVM is faster than DBZ by almost two orders of magnitude, and provide experimental results that validate the theoretical complexity analyses.

Creative telescoping algorithms compute linear differential equations satisfied by multiple integrals with parameters. We describe a precise and elementary algorithmic version of the Griffiths-Dwork method for the creative telescoping of rational functions. This leads to bounds on the order and degree of the coefficients of the differential equation, and to the first complexity result which is single exponential in the number of variables. One of the important features of the algorithm is that it does not need to compute certificates. The approach is vindicated by a prototype implementation.

We show that linear differential operators with polynomial coefficients over a field of characteristic zero can be multiplied in quasi-optimal time. This answers an open question raised by van der Hoeven.

We study tight bounds and fast algorithms for LCLMs of several linear differential operators with polynomial coefficients. We analyse the arithmetic complexity of existing algorithms for LCLMs, as well as the size of their outputs. We propose a new algorithm that recasts the LCLM computation in a linear algebra problem on a polynomial matrix. This algorithm yields sharp bounds on the coefficient degrees of the LCLM, improving by one order of magnitude the best bounds obtained using previous algorithms. The complexity of the new algorithm is almost optimal, in the sense that it nearly matches the arithmetic size of the output.

We provide algorithms computing power series solutions of a large class of differential or q-differential equations or systems. Their number of arithmetic operations grows linearly with the precision, up to logarithmic terms.

The long-term goal initiated in this work is to obtain fast algorithms and implementations for definite integration in Almkvist and Zeilberger's framework of (differential) creative telescoping. Our complexity-driven approach is to obtain tight degree bounds on the various expressions involved in the method. To make the problem more tractable, we restrict to bivariate rational functions. By considering this constrained class of inputs, we are able to blend the general method of creative telescoping with the well-known Hermite reduction. We then use our new method to compute diagonals of rational power series arising from combinatorics.

We address complexity issues for linear differential equations in characteristic p >;0: resolution and computation of the p-curvature. For these tasks, our main focus is on algorithms whose complexity behaves well with respect to p. We prove bounds linear in p on the degree of polynomial solutions and propose algorithms for testing the existence of polynomial solutions in sublinear time Õ(p1/2), and for determining a whole basis of the solution space in quasi-linear time Õ(p); the Õ notation indicates that we hide logarithmic factors. We show that for equations of arbitrary order, the p-curvature can be computed in subquadratic time Õ(p1.79), and that this can be improved to O(log(p)) for first order equations and to Õ(p) for classes of second order equations.

International audience

Efficient algorithms are known for many operations on truncated power series (multiplication, powering, exponential, ...). Composition is a more complex task. We isolate a large class of power series for which composition can be performed efficiently. We deduce fast algorithms for converting polynomials between various bases, including Euler, Bernoulli, Fibonacci, and the orthogonal Laguerre, Hermite, Jacobi, Krawtchouk, Meixner and Meixner-Pollaczek.

It is classical that univariate algebraic functions satisfy linear differential equations with polynomial coefficients. Linear recurrences follow for the coefficients of their power series expansions. We show that the linear differential equation of minimal order has coefficients whose degree is cubic in the degree of the function. We also show that there exists a linear differential equation of order linear in the degree whose coefficients are only of quadratic degree. Furthermore, we prove the existence of recurrences of order and degree close to optimal. We study the complexity of computing these differential equations and recurrences. We deduce a fast algorithm for the expansion of algebraic series.

Linear systems with structures such as Toeplitz-, Vandermonde-or Cauchy-likeness can be solved in O~(α2n) operations, where n is the matrix size, α is its displacement rank, and O~denotes the omission of logarithmic factors. We show that for Toeplitz-like and Vandermonde-like trices, this cost can be reduced to O~(αω--1 n), where ω is a feasible exponent for matrix multiplication over the base field. The best known estimate for ω is ω< 2.38, resulting in costs of order O~(α1.38n). We also present consequences for Hermite-Padé approximation and bivariate interpolation.

We study the complexity of computing one or several terms (not necessarily consecutive) in a recurrence with polynomial coefficients. As applications, we improve the best currently known upper bounds for factoring integers deterministically and for computing the Cartier–Manin operator of hyperelliptic curves.

We consider two kinds of problems: the computation of polynomial and rational solutions of linear recurrences with coefficients that are polynomials with integer coefficients; indefinite and definite summation of sequences that are hypergeometric over the rational numbers. The algorithms for these tasks all involve as an intermediate quantity an integer N (dispersion or root of an indicial polynomial) that is potentially exponential in the bit size of their input. Previous algorithms have a bit complexity that is at least quadratic in N. We revisit them and propose variants that exploit the structure of solutions and avoid expanding polynomials of degree N. We give two algorithms: a probabilistic one that detects the existence or absence of nonzero polynomial and rational solutions in O(√N log2 N) bit operations; a deterministic one that computes a compact representation of the solution in O(N log3 N) bit operations. Similar speedups are obtained in indefinite and definite hypergeometric summation. We describe the results of an implementation.

We investigate polynomial solutions of homogeneous linear differential equations with coefficients that are polynomials with integer coefficients. The problems we consider are the existence of nonzero polynomial solutions, the determination of the dimension of the vector space of polynomial solutions, the computation of a basis of this space. Previous algorithms have a bit complexity that is at least quadratic in the largest integer valuation N of formal Laurent series solutions at infinity, even for merely detecting the existence of nonzero polynomial solutions. We give a deterministic algorithm that computes a compact representation of a basis of polynomial solutions in O(Nlog3N) bit operations. We also give a probabilistic algorithm that computes the dimension of the space of polynomial solutions in O(√Nlog2N) bit operations. In general, the integer N is not polynomially bounded in the bit size of the input differential equation. We isolate a class of equations for which detecting nonzero polynomial solutions can be performed in polynomial complexity. We discuss implementation issues and possible extensions.

Many polynomial factorization algorithms rely on Hensel lifting and factor recombination. For bivariate polynomials we show that lifting the factors up to a precision linear in the total degree of the polynomial to be factored is sufficient to deduce the recombination by linear algebra, using trace recombination. Then, the total cost of the lifting and the recombination stage is subquadratic in the size of the dense representation of the input polynomial. Lifting is often the practical bottleneck of this method: we propose an algorithm based on a faster multi-moduli computation for univariate polynomials and show that it saves a constant factor compared to the classical multifactor lifting algorithm.

The transposition principle, also called Tellegen's principle, is a set of transformation rules for linear programs. Yet, though well known, it is not used systematically, and few practical implementations rely on it. In this article, we propose explicit transposed versions of polynomial multiplication and division but also new faster algorithms for multipoint evaluation, interpolation and their transposes. We report on their implementation in Shoup's NTL C++ library.