Creative telescoping using reductions

The technique of creative telescoping is a powerful tool for proving and finding combinatorial identities and the computation of parametric integrals that involve special functions. The name was coined by van der Poorten [77] in relation with Apéry's irrationality proof of

, but the first systematic algorithms were only developed about one decade later by Zeilberger and followers [90, 4, 91, 62, 85, 86]. Early precursors are [35, 89, 42].

Let us briefly recall the main idea behind creative telescoping on an extremely simple example: assume that we wish to prove the well known identity

Summing this relation from

until

, the right hand side becomes a telescoping sum, and we obtain

. Since

, the recurrence relation

proves

More generally, for the evaluation of a sum of the form

, the idea is to construct a relation of the form

where the left-hand side does not depend on

and the right-hand telescopes when summing over

. If we manage to construct such a relation, then the creative telescoping process provides us with a non-trivial linear recurrence relation for

. Writing

for the operator that sends

, the difference operator

that annihilates

in (1.1) is called a telescoper, and the function

its corresponding certificate.

A similar idea works for integrals of the form

instead of discrete sums. In that case, we need a relation of the form

instead of (1.1). Again,

is called the telescoper and

the certificate. Several other variants can be considered as well, such as

-difference equations, multiple sums and integrals, etc.

After Zeilberger's initial impetus, the development of creative telescoping focussed on making implementations more efficient, while guaranteeing termination and widening the class of functions on which the theory can be applied. Whereas early implementations were particularly successful for first or low order recurrence relations (e.g. hypergeometric summation), efficient algorithms soon appeared for more general holonomic functions and D-finite ideals in Ore algebras [26, 29, 27, 63]. For a historical perspective on these advances and further references, we refer to [28].

One major application of creative telescoping is proving functional identities. In the purely differential setting, it should be noticed that alternative methods have been developed for larger classes of non-linear equations [32, 33, 83, 75, 46, 49].

1.2.From Hermite reduction to reduction-based telescoping

Despite the above successes, the efficiency of creative telescoping remained a major issue. Part of the problem was due to the fact that the construction of relations of the form (1.1) or (1.2) often relied on brute force computations rather than a better insight into their precise shape.

The next major step was the introduction of reduction-based algorithms as a remedy to these problems. The concept of a reduction goes back to the works of Ostrogradsky and Hermite on the integration of rational functions [74, 44]. Hermite reduction is nowadays a fundamental tool in symbolic integration and generalizations for algebraic, hyperexponential and hypergeometric functions were introduced with this kind of application in mind [88, 30, 2]. It is interesting to notice that the importance of Hermite reduction for integrals that depend on parameters was understood quite early [76]; see also [39].

Hermite reduction was first related to creative telescoping in [41]; the complexity perspective appeared in [13]. Initially, only quite restricted types of functions could be treated this way and it became a major challenge to develop reduction-based algorithms for similarly general classes of functions as for the other approaches. This triggered a lot of activity in recent years [24, 14, 67, 22, 34, 23, 54, 68], culminating with an algorithm for Fuchsian differential equations [21]. This last result made it plausible that a fully general algorithm in the differential case might be quite complex, since it might involve elaborate arguments in order to deal with irregular singularities. Fortunately this concern turned out to be unwarranted: last year we constructed suitable reductions for arbitrary differential equations in a purely algebraic way [50]. Our technique has recently been further improved in [16]. In the present paper, we continue its development and extend it to mixed differential-difference equations.

Before we proceed, it is useful to recall how reductions show up in relation with creative telescoping. Assume that we wish to find a linear differential equation over

for the function

, where

. Hermite reduction in

provides us with a unique decomposition of any

as a sum

where

admits only simple poles in

and

. When restricting ourselves to functions

whose poles in

belong to a fixed finite set (i.e.

for some fixed square-free polynomial

), it follows that

is “confined” to a finite dimensional

-vector space

. Indeed, for some algebraic extension

that contains all roots of

, the reduction

is always an

-linear combination of the fractions

, where

runs over the set of roots of

. In our case, we simply take

to be the square-free part of the denominator of

Now the crucial observation is that

, etc. all belong to the above finite dimensional vector space, whence there exists a non-trivial relation

with

. By construction, there exists an

with

for

. It is not hard to check that the reduction

commutes with partial derivation with respect to

, whence (1.4) implies

This new relation is of the desired form (1.2). Another advantage of this method with respect to previous ones is that the computation of certificates is possible, but only optional. Since certificates tend to be more voluminous than the telescopers themselves, this is very interesting from a complexity point of view.

In order to generalize the above reduction-based approach to functions

that satisfy higher order differential equations, one essentially needs to define a confined reduction on a so called narrow

-submodule

that contains

. In the differential case, such narrow modules are still of the same form

. In the difference case, the set of allowed poles are translates of a finite number of points and their orders need to be bounded in a way that will be made precise.

1.3.On the choice of our setup

The theory of creative telescoping has become very wide, with many variants that cover specific cases more or less efficiently. For the exposition in this paper, we therefore had to make several choices concerning the setup. Let us briefly comment on these choices.

Non-commutative operator algebras. It was recognized by Zeilberger that the framework of holonomic functions is particularly suitable for the development of creative telescoping [90]. The reason is that there exists a systematic elimination theory for the equations satisfied by such functions that can be regarded as a non-commutative counterpart of the theory of Gröbner bases.

Elimination methods for differential equations actually predate Buchberger's famous algorithm [19]. One line of development is due to Riquier–Janet–Thomas [79, 56, 87] and another one to Ritt–Kolchin [80, 61]. The most “effective” variant of the latter theory [81, 17] essentially covers both commutative Gröbner basis theory, its skew differential counterpart for Weyl algebras, and non-linear generalizations that specialize to triangular system type solving in the algebraic case. Yet another line of development is due to Ore, who introduced a systematic theory for non-commutative polynomials [72, 73]. A more direct inspiration for Zeilberger's work were developments in the theory of D-modules in the seventies [9, 10, 12] and the development of non-commutative Gröbner basis theory in the computer algebra community in the eighties and later [40, 84, 60, 71, 64, 26, 70].

In this paper, we have chosen to work over differential-difference algebras generated by a finite number of coordinates

, derivations

, shift operators

, and

-difference operators

. All derivations

and automorphisms

are assumed to commute pairwise. Our setting is not the most general possible, but it has a natural geometric interpretation, as well as a non-linear counterpart [45]. The commutation of the derivations is not really essential, except for the main derivation that will be eliminated. It should also be quite straightforward to replace the algebra

in section 6.4 by more general algebras of solvable type [60, 64], Ore algebras [26], or G-algebras [70, 43]. On the other hand, when telescoping with respect to a difference operator it is essential that singularities are moved in a bijective way. It is therefore natural to restrict one's attention to homographies. Then, modulo a change of variables, homographies always operate like shift operators or

-difference operators: see section 5.1. Other interesting generalizations would concern nested extensions by solutions of linear differential-difference equations as in [78, 82].

Scalar equations versus first order systems. An eternal dilemma when dealing with linear differential or difference equations is whether one should privilege scalar higher order equations or first order systems. Thanks to the cyclic vector lemma [3], both approaches are essentially equivalent, so the choice is mainly a matter of taste. The best solution is probably to work with the representation that was used for the input, whenever possible. Our personal taste is slightly in favour of scalar equations, although we opted for first order systems in [50]. In [16], the authors preferred to work with scalar equations instead. In the present paper, we mainly consider scalar equations, but also outline how to adapt the theory to first order systems. At the end of the day, it seems that both approaches are more or less equally diligent.

Algebraic extensions versus rational counterparts. From the early days of symbolic computation, algorithms that require computations in algebraic extensions have acquired the reputation of being slow. Starting with Trager's work on the integration of rational functions [88], it has therefore become common practice to develop “rational” counterparts for reduction algorithms that avoid explicit computations in algebraic extensions. Nevertheless, from a complexity point of view, it is not so clear that working in algebraic extensions is necessarily that bad: computing in an extension of degree

is roughly

times more expensive, but typically deals with

conjugate roots in a single pass. Although the development of rational counterparts remains an interesting topic, we therefore do not regard such algorithms as intrinsically better. In this paper, we have chosen to systematically work in algebraic extensions whenever appropriate, which is somewhat simpler from a conceptual standpoint.

Characteristic zero. It will be convenient to assume that all fields considered in this paper are of characteristic zero, although most results can be generalized in a straightforward way to fields of (sufficiently large) positive characteristic.

1.4.Structure of the paper and outline of the main results

The main purpose of this paper is to present reduction-based algorithms for creative telescoping in the above setting of differential-shift-

-difference equations and to derive complexity bounds for several of these algorithms. In the purely differential case, our algorithms are completely systematic. In the presence of difference operators it is well known that indefinite summation and integration under the integral sign do not always preserve D-finiteness. Nevertheless, we provide a sufficient condition (explicit telescopability) for our algorithms to work, as well as a reduction-based algorithm that is able to check D-finiteness in a particular case. Let us briefly outline the structure of the paper and summarize our contributions.

Abstract reductions. For a gentle introduction to reduction-based creative telescoping, we refer to [34, Section 1.2]. In section 2, we start with a more abstract presentation of the basic theory that covers all cases that occurred in the literature so far. We also introduce the notion of “local reductions” and describe various related constructions on an abstract level.

The Lagrange identity. Let

be a field. One major new idea in [16] is the consideration of reductions with respect to general differential operators

instead of just the derivation as in (1.3). More precisely, a reduction with respect to

is a

-linear projection

such that

for all

. Now in the case of a rational function

, ordinary Hermite reduction provided us with a way to derive telescoping relations (1.5) from linear dependencies (1.4). In the case of a function

that satisfies a linear differential equation

, it is shown how to do something similar using reductions with respect to the adjoint operator

. The technical tool that makes this possible is the Lagrange identity [55, section 5.3]. We recall this formula in section 3 and generalize it to difference and matrix operators.

Differential reductions. In section 4, we show how to construct reductions for differential equations. This case was already covered in [50], but it is instructive to present it with the formalisms from sections 2 and 3 at hand. In sections 4.2 and 4.3, we start with a variant of the construction from [16], whereas section 4.5 explains the link with [50] and also clarifies how the Lagrange identity can be used to simplify the head/tail chopping process. We also show how new-style reductions with respect to the adjoint operator

give rise to old-style reductions in the sense of [50] and vice versa. In section 4.6, we finally consider mixtures of both reduction styles. This makes it possible to extend creative telescoping to the resolution of more general linear differential or difference equations that involve parameters: see Remark 6.7.

Difference reductions. The analogue construction of reductions for difference operators is the subject of section 5. Whereas derivations do not introduce new singularities, but only aggravate the order of poles, difference operators move singularities, but do not increase their severity. The idea behind the construction of differential reductions is to diminish the order of the poles of the function to be reduced. Likewise, difference reductions are based on the idea to shift singularities back until they are confined in a finite set. Another difference with section 4 is that the construction is somewhat less canonical due to the absence of a privileged section of orbits of singularities (see Remark 5.7), but this is of no consequence for the application to creative telescoping. Apart from that, the theory from section 4 naturally adapts to the difference setting.

Creative telescoping. After the formal introduction of the setting of DD-operator algebras in section 6.1, we next consider the application to creative telescoping. For suitable “telescopable” D-finite ideals

of a DD-operator algebra

, the theory from sections 4 and 5 allows us to define a computable confined reduction on a narrow submodule

. In section 6.4 we show how this can be used to construct telescopers. The algorithm is a straightforward generalization of the technique from section 1.2 and only relies on linear algebra. As an end-result, it produces a finite set of generators for a D-finite ideal of telescopers; in the favourable case when the reduction is “normal”, these generators actually form a Gröbner basis of the ideal of all telescopers. In the context of Gröbner basis, the algorithm is reminiscent of the FGML algorithm [37] and a variant of the one used in [16].

D-finiteness tests for telescoping ideals. In the purely differential case, it is a well-known fact from the theory of D-modules (reproved in section 7) that all D-finite ideals are telescopable, so the algorithm from section 6.4 systematically works. In the presence of difference operators, the telescopability of a D-finite ideals depends on the nature of its singularities. In section 7, we examine such singularities in more detail and derive sufficient conditions for the telescopability of a D-finite ideal with a given basis for the quotient module (a particular was treated before in [16, section 4.4]). Sometimes these “explicit telescopability” conditions are only met after a change of basis. For differential reductions, we present an algorithm to compute such a base change if it exists. Even when the explicit telescopability conditions do not hold, it sometimes still happens that the telescoping ideal is D-finite; in the last subsection 7.4, we present an algorithm that decides whether this is the case and that computes a Gröbner basis of the telescoping ideal if it is D-finite.

Complexity of rational function arithmetic. The last two sections are devoted to analyzing the complexity of reduction-based algorithms for creative telescoping. For this, we first recall basic complexity results for computations with rational functions in several variables, including partial fraction decomposition. When using the common dense representation, these complexity bounds are quite pessimistic, yet always polynomial if the number of variables is fixed. When allowing for probabilistic algorithms of Las Vegas type, it is more efficient to represent rational functions by the straight line programs (SLPs) [20], and fast dense univariate polynomial arithmetic can be used at every evaluation point.

Complexity of creative telescoping. In the last section we finally turn to the cost of creative telescoping itself. For this, we use a separate analysis for the complexity of the reduction process and the complexity of the FGLM-type algorithm for computing linear relations. The latter algorithm relies on standard linear algebra and its complexity analysis is straightforward in terms of the dimension of the vector space

of reduced functions. We next turn to the reduction process and restrict our attention to differential reductions. In this case, we prove a polynomial complexity bound, as well as a polynomial bound for

as a function of the input size. Altogether, this yields a polynomial time algorithm for creative telescoping, as promised in [50]. The complexity exponent is again extremely bad when using dense representations, but essentially cubic per evaluation point in the SLP model.

Remark 1.1. It is well known that it is easy to construct operators in

of small bit size that admit polynomial solutions with a large number of terms. For instance, the operator

admits

as a solution. For the design of polynomial time algorithms for operators in

it is therefore crucial to avoid the explicit computation of such solutions. For a similar reason, telescopers of minimal order (e.g. of the form

in (1.5) with minimal

) can admit very large coefficients. The telescopers computed by our polynomial time algorithm in section 9 are therefore not necessarily minimal. The above phenomenon also implies that the algorithms at the end of section 7 (for base changes to reach an explicitly telescopable form and creative telescoping for not explicitly telescopable ideals) do not run in polynomial time, in general.

1.5.Pending issues

Let us briefly list a few problems that have not been tackled in the present draft by lack of time. We intend to address them in a future version or in a follow-up paper.

Base changes for explicit telescopability. The base change algorithm to remove apparent singularities from the end of section 7.3 has been presented in the differential case only. We are still looking for existing work on this topic in the literature and a similar algorithm for the difference case.

Testing D-finiteness. The general algorithm for creative telescoping from section 7.4 (in absence of explicit telescopability) has also been presented in the differential case only. The difference analogue seems to require a bit more work, but it should be feasible to work out an algorithm along similar lines.

Fast arithmetic for rational functions. Several complexity bounds from section 8 are not very sharp or a bit sketchy. We intend to improve this section in the next version.

Complexity bounds for the difference case. The final complexity bounds in sections 9.2 and 9.3 have only been presented in the differential case (which still allows for difference operators among the

). If

and

in Proposition 7.7 do not grow too fast, then it should be possible to prove polynomial complexity bounds in the difference case along similar lines (this in particular covers the case when

for

). If

becomes large, then the defining equations of

tend to become large as well due to Proposition 7.2. The hope is that a polynomial complexity bound can somehow be derived from this.

Acknowledgments. We are grateful to Grégoire Lecerf for a few helpful thoughts concerning section 8.

2.Abstract reductions

2.1.Definitions and basic properties

Let

be a field and consider a linear map

on some

-vector space

. A reduction with respect to

is a linear map

that satisfies:

The condition R2 stipulates that

is a projection. It will be convenient to also introduce the complementary projection

involved in R1, called the remainder:

A normal reduction

with respect to

can be regarded as a projection of

onto a subvector space

that is isomorphic to the cokernel of

Remark 2.2. Sometimes, the linear map

is only defined on a

-vector space

that contains

. Our definitions easily extend to this case module the replacement of

in the conditions R1 and R4.

2.2.Local reductions

Let

be a projection and

. Given a reduction

with respect to a linear map

, we say that the reduction is local for the projection

Proof. Given

, there exists a

with

, whence

. Since

, we also have

. We conclude that

.□

Proof. Given

, there exists a

with

and

. Furthermore,

, whence

and

. We finally have

.□

2.3.Gluing local reductions together

Two projections

are said to be orthogonal if

. A family

of projections

is said to be orthogonal if its members are pairwise orthogonal. Such a family is said to be summable if

is finite for all

. In that case, the sum

with

is well defined for all

. Moreover,

is a projection that is orthogonal to each

and any

can be decomposed canonically as

Two local reductions

for

are said to be independent if

and

are orthogonal and for all

, we have

Such reductions necessarily commute: we have

, and

in a similar way. A family

of reductions for projections

is said to be independent if the

are pairwise independent and if

is summable.

Proof. For all

and

, we get from (2.1) that

, whence

. It follows that

, which proves (2.4). We also have

, whence

, and

. It follows that

. Since

, this shows that

is a reduction with respect to

. We also have

, so

is local for

Let us prove (2.5) by induction over

. To simplify notations, assume that

for all

, and write

for

. For

, we have nothing to prove. For

, the result follows from (2.1). If

, then let

and

. Notice that

for

. The induction hypothesis yields

, whence

2.4.Composition of reductions

Proof. Let us first assume that

. For each

, we have

. We also have

. The general case follows using an easy induction on

.□

Let

be two local reductions with respect to

for two orthogonal projections

. We have seen that Proposition 2.6 applies whenever (2.2) holds. If we only have

for all

, then the following still holds:

2.5.Normalization of confined reductions

Consider a confined reduction

with respect to

. Then

is a finite dimensional vector space and

a finite dimensional subvector space of

that we call the space of exceptional functions for

. Notice that

. Given any projection

with

and

, let

be defined by

Of course, we have

. Finally, for any

, the construction ensures that

whence

.□

This proposition also admits a local analogue. Assume that we are given a local reduction

with respect to

for the projection

. We say that

is locally confined if

is finite dimensional, and locally normal if

. Assume that

is locally confined and decompose

with

. Let

a linear mapping such that

is a projection of

(such a map is easily constructed from a basis for

). Then we define

Proof. Given

, let

. Then we have

. Letting

be such that

, we also have

, whence

, since

and

coincide on

. From

, we conclude that

, whence

is a reduction with respect to

The reduction is local for

, since

. By construction,

. For

and

, we finally have

.□

3.Lagrange's identity and generalizations

3.1.The original differential case

Let

be a differential field for the derivation

and let

denote the skew ring of linear differential operators over

. For instance, one may take

and

. Recall that the adjoint of a differential operator

is defined by

Given two indeterminates

and

, we write

for the set of

-linear combinations of products

with

. The following identity is due to Lagrange [66]:

Proof. Let us first prove the existence of

by induction on

. If

, then

, and we may take

. If

, then we write

with

of order

and

. By the induction hypothesis, there exists a

with

3.2.The difference case

Assume now that

is a difference field for the automorphism

and let

denote the corresponding ring of skew difference operators. For instance, one may take

and

with

. We recall that the adjoint of a difference operator

is defined by

Proof. We again prove the existence of

by induction on

. If

, then

, and we may take

. If

, then we write

with

of order

and

. By the induction hypothesis, there exists a

with

3.3.Matrix versions

Let

again be a differential field for

and consider a not necessarily commutative differential algebra

over

. One important example of such an algebra is the algebra of

matrices

. Proposition 3.1 naturally generalizes to the case when

admits coefficients in

, modulo a few precautions: in this setting, it is important to distinguish between usual operators

that operate on the left

where

are the coefficients of

. Similar precautions apply to the difference setting, and we have:

3.4.Twisting

Let

be a differential field for

and assume that

is the logarithmic derivative

of another element

. In fact, we may allow

to live in some abstract extension field

. Then we have

for all

and we call

the twist of

. More generally, for any

, we define the twist

. One may check that

and

for all

. Using twisting, it is possible to replace the operator

at the right hand side of (3.1) by any first order differential operator:

Proof. Let

be a formal solution of

in some extension field

. Let

and

. Then Proposition 3.1 provides us with

such that

This proposition again admits a matrix generalization: for

, it suffices to replace

and

by matrices

with

invertible and

. We also need to assume that

commutes with each coefficient

One interesting special case occurs as follows: let

, let

be monic of order

and assume that we wish to obtain a generalization of (3.1) with

instead of

in the right hand side. Then we lift

into an operator

by sending each coefficient

and take

This equation can be considered as the “desired” generalization of (3.1) with

instead of

in the right hand side.

As usual, Proposition 3.5 also admits a difference version. Let

be a differential field for

and let

be such that

. The twist of an operator

is again defined as

, so that

3.5.Another twisted variant

There are many more variants and generalizations of Lagrange's identity. Let us describe one more variant in the differential case. Let

be a possibly non-commutative

-algebra, where

is a field.

Proof. Writing

for the derivative of a differential operator

, and using the expansion

where

plays a similar role as

in (3.1), where

plays the role of

, and where the derivation

got replaced by a general monic operator

4.Differential reductions

4.1.Partial fraction decomposition and Hermite reduction

Let

be a field with algebraic closure

. Any rational function

admits a partial fraction decomposition

Given

, we will write

for the projection

and we call

the polar part of

. The number

stands for the order of the pole of

. We also write

for the projection

and call

the constant part of

The decomposition (4.1) has the advantage of being symmetric with respect to all points in

. Often, the polar part at infinity is combined with the constant part, in which case (4.1) becomes

We also denote

and write

for the projection

. Notice that

forms a ring for any

The Hermite reduction

with respect to the ordinary differentiation

with respect to

is defined by

It is not hard to check that this indeed defines a normal reduction on

with respect to

. Its restriction to

determines a normal reduction on

with respect to

, by Proposition 2.3.

4.2.Local reduction with respect to differential operators

Let us now consider an arbitrary non-zero linear differential operator

instead of

. For any

, it is well-known [55] that there exists an indicial polynomial

and a shift

with

Given

, let us now construct a local reduction

. Given

, let

and

be such that

In a similar way, we may construct a local reduction

at infinity. Given

, let

and

be such that

Proof. We adopt the notations from above. Whenever

, we let

and

. Whenever

, we also denote

and

By induction on

, let us first show that there exists a

with

and

; notice that this in particular implies ii, as well as

. If

, then we may take

. If

and

, then there exists a

with

and

, whence

. If

and

, then there exists a

with

and

. Hence

, so we may take

with

Again by induction on

, let us next show that

. If

, then we have

. If

and

, then

implies

. Consequently,

. If

and

, then

. This completes the proof that

is a reduction with respect to

Decomposing

with

, let us next show by induction on

that

. If

, then we have

and

. If

and

, then we may decompose

and

with

, after which

. If

and

, then

implies

. It follows that

. This shows i.

Now let

be such that

. The projections

and

are clearly orthogonal and it follows from ii that

and

; this shows iii. The last fact iv also follows from ii.□

Proof. Adopt the same notations as in the proof of Proposition 4.1. Given

, let us prove by induction on

that

. This is clear if

. If

and

, then the induction hypothesis yields

, whence

. If

and

, then the induction hypothesis directly yields

. This completes the proof that

. If

, then we also have

and

. Since

admits at most

roots, the cardinality bounds follow.□

Remark 4.3. Proposition 4.2 shows that

is locally confined for all

. This makes it possible to define a locally normal reduction

using Proposition 2.10. In order to make this construction fully effective, we still have to determine the space

. Let

We claim that

. Assume for contradiction that there exists a function

with

and chose

such that

is minimal with this property. Let

as in the proof of Proposition 4.1. Since we must have

, it follows that

. Therefore

provides a counterexample with

, contracting our minimality hypothesis.

Having proved our claim, one may use linear algebra to determine elements

with

and such that

form a basis for

. Setting

and

, we define

to be the “distinguished” projection of

such that

for all

. In order to apply Proposition 2.10, we finally construct the mapping

by sending each

and then extend

by setting

Since

can be arbitrarily large, we notice that the determination of

can be very expensive. In particular, there exists no polynomial time algorithm to compute them. We also notice that the construction of

is symmetric in the sense that

for all automorphisms

over

4.3.Global reduction with respect to differential operators

From Propositions 4.1 and 2.6, it follows that we can glue the local reductions at the finite points

together into a reduction

, defined by

Notice that Proposition 2.3 implies that

whenever

. Using Propositions 4.1-iv and 2.8, we may next define the global reduction

with respect to

Remark 4.4. Applying the above construction to the local reductions

from Remark 4.3 instead of

, we claim that we obtain a normal reduction. Indeed, given

, we get

for all

, whence

and

. Since

is locally normal, we conclude that

. We notice that

is again symmetric in the sense that

. It is also easily checked that

is locally confined for each projection

It is interesting to consider the restriction of the global reduction

to certain differential subrings of

and

. More precisely, given a subset

of poles, the set

forms a differential

-subalgebra of

. From Proposition 4.1-ii, it follows that

is stable under reduction. If

is the set of zeros of a monic separable polynomial

, then we observe that

coincides with

. In that case, Proposition 2.3 implies that

is stable under reduction as well. We notice that

is finitely generated by

as a

-submodule of

; for this reason we call it a narrow submodule of

Proof. Given

and

with

, we notice that

, whence

. Given

with

, it follows that

. In other words,

is stable under reduction and so is

, since the reduction commutes with all automorphisms of

over

Let us now turn to the dimension bound. Since

-linear, it suffices to show that

. Let

be defined as in Proposition 4.2 for all

. Then the proposition implies that

For any

, we have

by Proposition 4.1-ii, whence

. Since

by Proposition 4.2, it follows that

, whence

Remark 4.6. A bit of flexibility is possible regarding the choice of

in the definitions of the local reductions

. Taking

for

as in [16] may lead to slightly more confined reductions in Theorem 4.5. Conversely, by taking

somewhat larger, one may ensure that

for all

in Proposition 4.1-ii, so that

and

commute for all

in Proposition 4.1-iii. This makes it possible to treat the singularity at infinity in a more symmetric way and replace definition (4.6) by

Remark 4.7. Since the reduction

form Remark 4.4 is both normal and locally confined for each projection

with

, its restriction to

yields a normal confined reduction.

4.4.Back to reductions with respect to

Assume still that

is a non-zero differential operator. Let

be its order, let

be its dominant coefficient, and let

be an arbitrary monic separable polynomial such that

is contained in

Let us show how to construct a confined reduction

with respect to

. We start with the confined reduction

with respect

provided by Theorem 4.5. Now consider

By construction, we have

, whence there exists a

such that

. Using Proposition 3.1, it follows that

4.5.First order systems

Let us now turn to first order systems and study analogues for the theory and results from the previous subsections. This means that the operator

from section 4.4 is now a first order operator in

and we need to adapt the theory of sections 4.2 and 4.3 to the adjoint operator

; recall that

operates at the right on matrices in

. By Proposition 2.1, we also recall that the construction of a reduction

with respect to

is equivalent to the construction of a reduction

with respect to

, for any invertible matrix

. In particular, for the construction of a local reduction at

, it suffices to consider operators of the form

with

As the analogue for the operator

from sections 4.2 and 4.3, we will therefore start with an operator

and show how to construct the local reduction

with respect to

. It will be convenient to first put

into an even simpler form through multiplication by a suitable invertible matrix

. Let

. In analogy with [50], we say that an invertible matrix

is a tail chopper for

if we either have

and

, or

, and for some

Proof. At every iteration of this loop, we notice that

increases by one (as long as

) and it is easily checked that T1 and

are loop invariants.□

Remark 4.10. The above algorithm for the computation of tail choppers is somewhat reminiscent of Abramov's EG-elimination method [1].

Given a tail chopper

as above, let

be the denominator of

(so that

), and consider the operator

. We have

We consider

as the “preconditioned” version of

with respect to which we will now construct the local reduction

. From (4.9) it follows that the shift of

is simply

and its indicial polynomial is given by

By construction,

is invertible as a matrix in

, whence there exist at most

values of

where

. We are now in a position to formulate a counterpart of (4.3).

Although this definition works, it is slightly suboptimal in the sense that values of

for which

are “skipped” altogether even though the system

might actually be solvable. Given an arbitrary scalar matrix

, let

denote a pseudo-inverse such that

for any

. Also let

denote a matrix with

for all

. We may now replace (4.10) by

This “more confined” definition ensures that

instead of

. The local reduction

at infinity can be defined in a similar fashion, via the change of variables

Given an arbitrary first order operator

, this completes our construction of the local reductions

and

with respect to

. Proposition 4.1 and the theory from section 4.3 can be adapted mutatis mutandis. In particular,

and

are indeed local reductions for

and

that can be glued together into a global reduction

with respect to

. Assuming that the mappings

and

were taken to commute with all automorphisms of

over

, we again have

for all

. The analogue of Theorem 4.5 is:

As to the counterpart for section 4.4, assume that

with

. Let

be an arbitrary monic separable polynomial such that

, where

. Given a formal column vector

with entries

, the set

admits the natural structure of a

-module for the derivation

with

. Given a row vector

, we define the reduction on

where

is the confined reduction on

with respect to

as constructed above. By construction, there exists a

with

, whence

Remark 4.12. Starting from a reduction

for

with

, one may also define a reduction

for

: given

, there exists a

with

, and we take

. Then there exists a

with

. Using Lagrange's identity the other way around, we obtain

, whence

. This shows that

is indeed a reduction with respect to

4.6.Reductions on D-modules with respect to general operators

We have seen how to construct reductions on

with respect to general operators and how to construct reductions on D-modules with respect to

. It is natural to ask whether we can construct reductions on D-modules with respect to more general operators.

Let us first show how to do this in the setting of first order systems from the previous subsection. More precisely, assume that

with

and

are as in the previous subsection. Given another matrix

, our aim is to construct a reduction on

with respect to

. The idea is to use an explicit variant of Proposition 3.5. Given a matrix

, we notice that

Now we may reinterpret

as a first order differential operator of dimension

, for which we may construct a confined reduction

using the theory from the previous subsection. Given

, this means that there exists a

with

. Setting

The above construction admits several variants. Let us briefly sketch what happens if the first order matrix operator

is replaced by a monic operator

. This time, we rather rely on Proposition 3.7 and more precisely on formula (3.5). Given a row matrix

and using the fact that

, we have

The idea is now to construct a confined reduction

with respect to

. This can be done by generalizing the theory of sections 4.2-4.4 to allow for matrix coefficients, along the lines of section 4.5. We next define the confined reduction on

as usual by

. Given

, there then exists a

with

, whence

5.Difference reductions

5.1.The field

as a difference field

Let

be a field and consider a birational map

from the projective line over

into itself. Such an map

is called a homography and is necessarily of the form

Given a homography

, the field

becomes a difference field for the automorphism

that postcomposes with

for all

, and

. More generally, if

was already equipped with an automorphism

, then we may extend

into an automorphism

by taking

where

acts coefficientwise on rational functions

. Inversely, we must have

for automorphisms

of this kind, and we call

the homography associated to

Two important special cases are the shift operators

and the

-difference operator

In fact, we claim that the general case can essentially be reduced to one of these two special cases via a change of variables

. Indeed, such a change of variables transforms difference equations in

and

into difference equations in

and

. Assuming that we are allowed to extend

with the roots of the characteristic polynomial of

, we may first put

in Jordan normal form and then normalize it through division by

. After this, we have

for some

For what follows, we will always assume that

has infinite order. In the case of a shift operator

, this means that we should have

; for a

-difference operator

, the number

should not be a root of unity.

5.2.Difference modules

One particularity of automorphisms

as above with respect to the ordinary differentiation

is that the sets

are no longer stable under

. On the other hand, the set

is stable under

, where

and

. It will be convenient to expand somewhat more on this observation and introduce a few more notations in addition to those from section 4.1.

The definitions extend to the case when

by taking

. Given

and a map

, we also denote

Each of these vector spaces also come with projections similar to the

and

that were already defined before. For example,

stands for truncation at order

. For

, it is also convenient to introduce

Given

, we write

if there exists an

such that

. This clearly defines an equivalence relation. We also write

if there exists an

such that

. Notice that

is a partial ordering on

by our assumption that

has infinite order. Given

, we introduce the subsets

and

For

, we understand that

, etc. We notice that

is a

-module, whereas

is a

-module. Since

is a fixed point for

, then we notice that

is even a

-module.

More generally, let

be an increasing mapping in the sense that

for all

. Then

is a

-module. If

, then the set

is actually a

-module. We say that

is symmetric if it is stable under all automorphisms of

over

. Similarly,

is symmetric if

is symmetric and

commutes with all automorphisms of

over

. In that case, the intersection

is a

-module and

-module, whenever

Assume that

and consider a

-submodule

with

. It is not hard to check that such a module is necessarily of the form

for some increasing

. We say that

is narrow if it is a finitely generated as a

-module. More generally, a finitely generated

-module of

of the form

is said to be narrow. Such a narrow submodule of

is always the intersection of

with a narrow submodule

for some symmetric

Consider a narrow submodule

that is finitely generated by

. Without loss of generality, we may replace each

by the finite set

and remove any elements

such that

for some

. This means that we may assume that

is finite and

for some function

. Now consider the map

defined by

5.3.Local and semi-local reduction with respect to shift operators

Let us now move our attention to a shift operator

, where

with

. We assume that

with both

and

. We took an operator with respect to

rather than

since, for the applications that we have in mind (see Proposition 3.2), the operator

will actually be the adjoint of an operator in

Since

, the action of

at infinity can be approximated to any order by a differential operator in

. In particular, there still exist an indicial polynomial

and a shift

(in fact,

, since

admits coefficients in

) such that

Consequently, we may define the local reduction

at infinity with respect to

in a similar way as in subsection 4.2.

Our next aim is to define a “semi-local” reduction

. Since

admits a singularity at

if and only if

admits a singularity at

for all

, we also need to take into account the behaviour at each of these shifted singularities. For this reason, the reduction

will only be “semi-local”.

Given

, we define the

-span of

to be the largest integer

such that

. If no such integer exists, then we set

. Assuming that

, consider

Proof. If

, then we have already noticed that

is a reduction with respect to

such that i and ii are satisfied. If

, then the proof follows the same scheme as in Proposition 4.1; for the sake of completeness, we give it in full. Throughout the proof, we denote

and

By induction on

, let us first show that there exists a

with

; notice that this in particular implies ii and also

. If

, then we may take

. If

and

, then there exists a

with

, whence

. If

and

, then there exists a

with

, where

. Hence

, so we may take

Again by induction on

, let us next show that

. If

, we have

. If

and

, then

implies

. It follows that

. If

and

, then we have

. This completes the proof that

is a reduction with respect to

Decomposing

with

, let us next show by induction on

that

. If

, then we have

and

. If

and

, then we may decompose

and

with

, after which

. If

and

, then

implies

. It follows that

. This shows i.

Now let

be such that

. Then the projections

and

are clearly orthogonal and it follows from ii that

and

; this shows iii. The next fact iv also follows from ii.□

Proof. Assume the notations from the previous proposition. Let us prove i by induction on

. If

, then

. If

and

, then

; since the induction hypothesis implies

, we obtain

. If

and

, then the induction hypothesis yields

Let us next prove ii, again by induction over

. If

, then we may take

. If

and

, then let

be such that

, and

. Then we have

, and

. If

and

, then let

be such that

, and

. We have

. Since

, we have

, and

. It thus suffices to take

5.4.Aligned downward reduction with respect to shift operators

Let us now turn to the problem of combining the semi-local reductions

where

for some subset

. Given a general subset

, the set

admits the partitioning

In this subsection we start with the case when

is aligned in the sense that

with

. In that case, we notice that

Although the construction from the previous subsection provides us with a local reduction

for the projection

, this reduction is sometimes suboptimal due to the fact that it may take a long time to reduce functions for which

is small. For this reason, we introduce an alternative semi-local reduction

We call

the downward reduction with respect to

for

. It is “less confined” than

but relies on the reduction

in order to reduce the polar parts for the singularities in

Let us first show that

is a reduction with respect to

. We clearly have

. By Propositions 5.1-ii and 5.2-iii, we notice that

for

, whence

and

. Then

equals

The reduction

is also local for

since

. The relation

shows that

. The properties iii and iv are shown in a similar way as for Proposition 5.1.□

Proof. With the same notations as in the previous proof, Proposition 5.2-i implies

for

, whence

. This proves i.

Proposition 5.2-ii implies the existence of

with

, and

, for

. In particular, this yields

. Now consider

with

. For

, we conclude that

As to iii, we have

, thanks to Proposition 5.2-iii. Proposition 5.2-iv similarly implies iv, since

for

.□

Proof. Let

and decompose it as

as above. From Proposition 5.3-ii, it follows that

. Notice that

if and only if

for

. Assuming that that

, it follows from Proposition 5.4-iv that

for

, whence

.□

Remark 5.6. In analogy with Remark 4.3, one may define a local normalization

. This time, it follows from Proposition 5.4 that the space

is spanned by

. Again the computation of a basis of

can be expensive and there exists no general polynomial time algorithm for doing so.

5.5.General downward reduction with respect to shift operators

Let us now assume that

is such that

is finite for every

. The equivalence relation

restricted to

then leads to a natural partitioning

By Propositions 5.3-i, 5.3-iii and 2.6, we get that

is a local reduction for

with respect to

. From Propositions 5.3-iv and 2.8, we deduce that

is a reduction with respect to

. We call

and

the downward reductions with respect to

for

and

. If

is symmetric, then

and

, by Proposition 2.3.

Remark 5.7. When taking

to contain a section of

with respect to the equivalence relation

, we ensure that

. Nevertheless, each

can only contain “half” of the set

, since

is assumed to be finite. This means that

cannot cover

, so

cannot truly qualify as a “global reduction”. In order to obtain a global reduction, one needs to combine downward and “upward” reductions (constructed from downward reductions with respect to the operator

). Even then, the obtained “global” reduction crucially depends on the choice of

. In fact, this kind of reductions will not be useful for our application to creative telescoping; we rather need reductions on suitable narrow submodules.

In order to obtain a counterpart for Theorem 4.5, let us now turn our attention to restrictions of downward reductions to suitable narrow submodules of

Proof. Let

be partitioned as in (5.3). Let us first show that the narrow

-submodule

is stable under downward reduction

. The module

is clearly stable under the projection

for any subset

. Since

is local for

, we also have stability under the reduction

. Now given

and

, we have

with

and

by Corollary 5.5. For any

, it follows that

and

. Since

commutes with all automorphisms of

over

, we also get the stability of

under

Given

, let us next prove the existence of some

with

. For each

, setting

, Proposition 5.4-ii yields an element

with

and

for all

, whence

. At infinity, there also exists a

with

, and

. Taking

, it follows that

and

admits coefficients in an extension field

of degree

, then we may replace

, since

and

for every automorphism

over

. Notice that this shows in particular that

is a reduction on

which also implies the bound on

. For each

, let

be as in (5.2) and notice that

satisfies

and

. Now given

, Proposition 5.3 implies

for all

, where

. It follows that

. We also have

, whence

. Since

, we conclude that (5.4) indeed holds.□

Remark 5.9. When using the locally normal reductions

from Remark 5.6 instead of the reductions

and similarly for the local reduction at infinity, it can be checked that the resulting reduction

in Theorem 5.8 is both confined and normal.

5.6.Back to reductions with respect to

Assume now that

is a difference operator of order

with

. Consider a narrow

-submodule

, where

is a symmetric finite set,

a symmetric function, and

. Assume that for every zero

of multiplicity

, we have

admits the natural structure of a difference module over

for the shift operator

with

and it suffices to check that

. Now for any root

of multiplicity

, we have

, whence

, thanks to our assumption (5.5). We conclude that

, as desired. Notice that

is a formal solution to the equation

Let us now show how to construct a confined reduction

with respect to

. Let

denote the restriction of the directed reduction at

with respect to

as constructed in the previous subsection. Now consider

By construction, we have

, whence there exists a

such that

. By Theorem 5.8, we actually have

. Using Proposition 3.2, it follows that

Remark 5.11. The smallest set

and map

that satisfy our hypotheses are given by

and

. Assume that

and

are taken this way. If

for all

, then it follows that

for all

, whence

. However, in the other extreme case when

is aligned, we only have

for all

and the growth of

can be quadratic in

. This larger growth is due to the fact that we required functions of small span to be reducible fast in section 5.4. If we drop this requirement, then linear growth can be restored by changing the definition (5.1) to

where

and

. Compromises between both definitions are also possible, which should ultimately make it possible to systematically gain a factor

with respect to some of the complexity bounds later in this paper. However, the details are technical, so we reserve them for a future work.

5.7.The case of

-difference operators

The theory of the previous section naturally adapts to

-difference operators

(where

is not a root of unity), except that both zero and infinity now need to be treated apart.

So let us start with the construction of the local reductions

and

for

. This time, the behaviour of

at zero and at infinity is given by

where

and

(in fact, we have

and

). Since

is not a root of unity, we notice that

admits at most one solution for any

The local reductions

at other points

are constructed in a similar way as in subsection 5.3. It is not hard to check that the

with

are indeed local reductions that satisfy similar properties as in Proposition 5.1. Given a subset

, this allows us to glue them together into a reduction

in a similar way as in subsections 5.4 and 5.5. Mutatis mutandis, this leads to the following analogue of Theorem 5.8:

The results from subsection 5.6 also naturally adapt to the

-difference setting.

5.8.First order systems

Let us now outline how to adapt the theory of this section to first order systems. We assume that

for

; the case of

-difference operators can be treated in a similar way. Concerning the analogue of section 5.3, let us start with an operator

, where

is a monic separable polynomial and where

is an invertible matrix. Recall that

operates at the right on matrices in

The local reduction at infinity

is constructed in a similar way as for differential operators. The semi-local reductions

at other points

are somewhat easier to construct in the sense that there is no need for any auxiliary tail choppers. Given

, we define

by induction over

As in section 4.5, this definition can be further optimized through the use of pseudo-inverses:

The remainder of the theory and result from sections 5.3, 5.4, 5.5, and 5.6 can be adapted mutatis mutandis. In particular, for every symmetric finite set

, we obtain a reduction

that satisfies the following analogue of Theorem 5.8:

6.Creative telescoping

6.1.D-finite ideals of DD-operator algebras

Let

be a field that is equipped with

pairwise commuting operators

, together with derivations

and automorphisms

such that for each

we have either

and

. Given such a field

, we will say that the skew polynomial ring

is a DD-operator algebra. In practice, we usually have

for some constant field

, and each

is given by either one of the formulas

Modulo a homographic change of variables (and quadratic algebraic extensions), we have seen in section 5.1 that we may further reduce to the case when each

is either a shift or a

-difference operator. DD-operator algebras of this special type will be called standard.

A (left) ideal

is said to be D-finite if

is finite dimensional as a vector space over

. Given an

-module

and a “function”

, we say that

is D-finite if its annihilating ideal

is D-finite. We say that an ideal

is reflexive if for every

and

, we have

It is well known [40, 60, 71, 64, 26, 70] that the theory of Gröbner bases generalizes to skew polynomial rings such as

. As usual, this assumes a total ordering

on the set of “monomials”

that refines the divisibility relation on

. Let

denote this divisibility relation. Given a finite subset of

, a Gröbner basis

for the left ideal generated by this subset can be computed using a non-commutative version of Buchberger's algorithm. The initial segment

of reduced monomials for

under the Gröbner stairs for

admits the usual property that

forms a basis for the vector space

. Given

, the reduction of

with respect to

yields

with

6.2.Matrix representations of operators

Let

be a reflexive D-finite ideal of

. Consider any basis

and let

be the column vector with entries

(it will be convenient to call such a column vector

a basis of

as well). Given an operator

, the action of

gives rise to a matrix

with

are such that

, then we may compute the

-th column of

by reducing

with respect to

and writing the result as a linear combination of

Proof. Assume for contradiction that there exists a non-zero row vector

with

. Then we would get

. Writing

for the column matrix with entries

, we thus obtain

and

, which yields a non-trivial relation between the

.□

Since both derivations commute and the entries of

form a basis of

, this yields the first relation. The two other relations are proved in a similar way.□

Given

, we say that

is a cyclic vector for

form a basis for

. If

and

contains an element

with

, then it is well known that such a cyclic vector

always exists and can be computed [3, 25, 15]. These results can naturally be adapted to the case when

is a shift operator or a

-difference operator such that

contains an element

such that

is infinite. We will say that

acts non trivially on

if we are in either one of these situations. In the terminology of [6, section 5.3], this implies that

is simple, whence [6, Corollary 5.3.6] yields a way to compute cyclic vectors. Given a cyclic vector

, the left ideal of

of operators that annihilate

is principal and its unique monic generator has order

. In terms of matrices, these properties may be restated as follows:

6.3.Narrowness

where

is the constant field of

. Moreover, if

, then we require that

is either a shift operator

with

or a

-difference operator

, where

is not a root of unity. As usual, we set

, where we understand that

. We also write

and

Recall that a

-submodule

is said to be narrow in

if it is stable under

and finitely generated as a

-module. Given such a module

, we have shown in section 5.2 that there exists a finite symmetric set

and a symmetric

such that

. Let

be a fixed basis of

. A

-submodule

is said to be narrow in

if it is of the form

for a narrow submodule

Given a narrow submodule

as above and

, we have seen in sections 4.5 and 5.8 how to construct a computable confined reduction

with respect to the first order operator

as well as a corresponding computable confined reduction

with respect to

. Alternatively, one may rely on the constructions of reductions with respect to scalar linear differential and difference operators

, as we will detail now.

Let

be the entries of the basis

. We say that

is cyclic if

for

. Such a basis always exists by Proposition 6.3 and it comes with a unique monic operator

of order

with

. Multiplying

with the lcm of the denominators of its coefficients, we obtain a new operator

with

and

Assume that

and that we are given a narrow submodule

as above with respect to a cyclic basis

. Then we have shown in sections 4.2, 4.3 and 4.4 how to construct a computable confined reduction

with respect to

as well as a corresponding computable confined reduction

with respect to

Assume now that

and that we are given a narrow submodule

as above with respect to a cyclic basis

. Assume also that

for all

, where

denotes the set of fixed points of

. Then we have shown in section 5 how to construct a computable confined reduction

with respect to

as well as a corresponding computable confined reduction

with respect to

, where

. For

, we notice that

for some

-vector space

of dimension at most

. Given

with

and

, we may write

with

and

for

. By induction on

, we may then define a confined reduction on

. We have

Assume that there exists a

-submodule

that comes with a (computable) confined reduction

with respect to

. Then we say that

is telescopable in

for the chosen basis

, and we call

the associated reduction. By what precedes, this is the case whenever there exists a

-submodule

that is narrow in

6.4.Creative telescoping

With the notations from the previous subsection, given a function

in some

-module, the set

of telescopers for

with respect to

is defined by

From this definition it is immediately apparent that

forms an ideal of

. Given a telescoper

, an element

such that

is called a certificate for

Assume now that

is telescopable in

, with associated reduction

. We say that a function

is primitive if for each

, there exists an

with

. We write

for the set of such functions. Given

and an operator

with

, we notice that

implies the existence of some

with

, whence

and

. In particular, given an operator

with

and

, we have

Since

lives in the finite dimensional

-vector space

, this yields a way to determine telescopers through linear algebra, simply by searching

-linear relations between the

where

runs over

. For instance, taking

, one obtains a telescoper in

, which generalizes the introductory example (1.4). Doing this for each

shows in particular that

is D-finite.

In our DD-algebra setting, it is also natural to perform the linear algebra incrementally following the term order

, as in the FGLM algorithm [37]. This leads to the following algorithm to compute a D-finite ideal contained in

Algorithm 6.1

Input: a telescopable D-finite ideal with associated reduction and

Output: the reduced Gröbner basis for a D-finite ideal included in

monomials remaining to be treated

Gröbner basis being constructed

monomials under stairs

while do

Let be minimal for and set

if is not reducible with respect to then

Compute

(If , then we may also take )

if , then

else

return the reduced Gröbner basis for the ideal generated by

Proof. Throughout the algorithm, we observe that the elements

with

are

-linearly independent. Consequently,

never exceeds

, which ensures termination of the algorithm. By (6.3), it is also clear that we only insert telescopers in

. Let us now prove the bound (6.4).

Throughout the algorithm, it is easily verified that

contains only monomials that are smaller than

. It follows that the leading monomial of

, when inserting a new element into

. This in turn implies that, at the start of each iteration of the main loop, the set

contains precisely those monomials below

that cannot be reduced with respect to

. At the end of the loop, this means that

consists exactly of those monomials that cannot be reduced with respect to

. Consequently,

Assume finally that

is a normal reduction. Then it is straightforward to verify that (6.3) becomes an equivalence. In particular, given a telescoper of the form

, we get that

, whence

for all

, since the elements

with

are linearly independent. In other words,

, whence

. At the end of the main loop, this also means that the S-polynomials of any two elements in

necessarily reduce to zero; in other words,

is already a reduced Gröbner basis, so we may directly return it after the main loop.□

Remark 6.5. The important property of

is that

. If

, then it follows that we may also take

instead of

, since

Remark 6.6. Reduction-based algorithms for creative telescoping do not require the computation of certificates. Nevertheless, it is not hard to modify the algorithms such that certificates are computed along with the telescopers themselves.

Remark 6.7. It is possible to generalize the theory of this section to the case when

is replaced by a more general operator

that commutes with

. In the differential case, it essentially suffices to replace the confined reduction

with respect to

by a confined reduction with respect to

, as constructed in section 4.6. In the difference case, such reductions can be constructed along similar lines.

7.D-finiteness tests for telescoping ideals

7.1.Location of singularities

Adopt the same notations as in section 6.3. Let

be a fixed basis of

with entries

. Given an operator

, we define its order

and with respect to the basis

Proof. Assuming for contradiction that

and

, we have

, whence Proposition 6.2 implies

whence

. Assuming that

, it follows by induction on

, we get that

implies

. Since

is finite, we must have

for some

. If

, it follows similarly that

implies

, whence the existence of some

with

7.2.Orders at singularities

Let

be the monoid generated by the

with

and

. Similarly, let

denote the multiplicative monoid generated by the

such that

Proof. Let

be sufficiently large such that

for all

with

. This implies that

. Let

. Since

and

is invertible, we get

Since

for

, we have

in the above formula. By our choice of

, we also have

, whence

Proof. Let

be sufficiently large such that

. Since

and

is invertible, we have

7.3.Explicit telescopability

Proof. Assume that

is finite and let

, where

. Notice that

is stable under both

and

. Given

, Proposition 7.1 implies that

. Given

, this means that

Example 7.6. If

is a standard DD-operator algebra with

, then

is finite if and only if

does not depend on

for every

with

. In particular, we may take

in the proof.

is finite, with

, and the restrictions

and

satisfy

. It follows that

is a narrow

-module of

. In particular, the

-module

comes with a reduction

with respect to

We claim that

is stable under

. The stability under

is clear. Let

and

with

, so that

. Then for any

, we have

Given an iterated derivative

and

, we next observe that there exist integers

with

By Proposition 7.3, there are only a finite number of points

where

and

at these points. In other words, setting

, there exists a finite dimensional

-vector space

with

. We extend

to a confined reduction

with respect to

by setting

for all

and

. Notice that

is still stable under

, since

commutes with

.□

Example 7.8. Consider the case when

is a standard DD-operator algebra with

. If

and

, then it follows that the denominator of

as a rational function in

should be a power product of linear forms

with

and

. If this is not the case, then one may try replacing the operator

by its inverse

, for which one might have more luck. If either

is a

-difference operator, and the other one is a shift operator, then the denominator of

should be a polynomial in

. If both

and

are

-difference operators, then the denominator of

should be a power product of terms

with

, and

We say that

is explicitly telescopable in

with respect to the basis

if it satisfies the sufficient conditions from Propositions 7.5 or 7.7. These conditions present the advantage that they are formulated exclusively in terms of the singularities of

. For a fixed basis

, they are also easy to check. It remains to investigate how to compute a basis with respect to which

is explicitly telescopable, whenever such a basis exists.

We will restrict our attention to the case when

; an algorithm can probably be developed for the difference case as well, but we will leave this for a future work. Given a singularity

, let us show how to compute a basis

with

whenever such a basis exists. We outline the algorithm over

, but it is not hard to adapt it to work over

. The computed basis is no longer cyclic; with additional effort, it can be shown that there also exist cyclic bases

with

We start with the computation of a basis of formal local transseries solutions of the system

, given in the form of a fundamental matrix

. If there exist local solutions that are not Laurent series, then the singularity

clearly cannot be removable. Otherwise, the truncation

at a sufficiently high order

gives rise to a basis

such that

admits no singularity at

Unfortunately,

may introduce new singularities, so this algorithm needs to be tweaked a little more. We replace

with a new approximation at a order

that can be written as a product of diagonal matrices with entries of the form

with

, invertible constant matrices in

, and upper/lower triangular matrices of the form

. It is not hard to compute such an approximation. By construction, its inverse admits a factorization of the same form, whence it introduces no new singularities.

Even small equations such as

may admit solutions with a high valuation

. It follows that the above algorithm for changing bases may be very expensive and its cost is generally not polynomial in the input size.

7.4.Testing D-finiteness

is explicitly telescopable in

with respect to the chosen basis

and associated reduction

, then the telescoping ideal

is certainly D-finite for any

. Moreover, Algorithm 6.1 allows us to efficiently compute a Gröbner basis for a D-finite subideal of

; with more computational effort, one may even compute a Gröbner basis for

itself.

is not explicitly telescopable in

, even after a suitable change of basis, then telescoping ideals

are generally not D-finite, although they exceptionally might be. In this subsection, we will show how to modify Algorithm 6.1 in order to test D-finiteness of

and to compute a Gröbner basis for

if so. We will restrict our attention to the case when

, but a similar approach is likely to work in the difference case as well. It will also be convenient to assume that for any

and

with

, we have

for all

. This is the case for standard DD-operator algebras.

Assume that we have fixed a cyclic basis

with entries

and let

be the vanishing operator of

. Denote

, and

. Then

is a

-submodule of

that contains

, but that is not necessarily narrow in

Let

be the normal global reduction on

with respect to

as constructed in Remark 4.7. The construction from section 4.4 still works and leads to a corresponding normal reduction

. The restrictions

and

of these reductions to

and

are also normal. We call

the associated normal reduction to

(and for the chosen basis

). For any subset

, we notice that

. In particular,

for all

Proof. Let

be the row vector with entries

. Then

implies

. The result now follows from (4.8).□

Proof. Let

be such that

and let

be such that

and

. Modulo the replacement of

, we may assume without loss of generality that

for all

. Recall that

are pairwise distinct.

Consider the sequence

with

and

for all

. Let us show by induction on

that

and

for all

. These assertions hold by assumption for

, so assume that

. By Lemma 7.9, we have

For each

, the assumption that

therefore implies

. Assume for contradiction that

. Let

be such that

Since

and

, it follows that

. From

, and the normality of

, we also get that

. We conclude that

, which contradicts the induction hypothesis.

We claim that

are

-linearly independent. Indeed, given a non-trivial relation

with

, we would obtain the contradiction

We are now in a position to state the adaptation of Algorithm 6.1 to the case when the ideal

is not necessarily telescopable.

Algorithm 7.1

Input: a D-finite ideal with associated normal reduction and

Output: the reduced Gröbner basis for if is D-finite and otherwise

, ,

while do

Let be minimal for and set

if is not reducible with respect to then

if then let

else decompose and let

Let be such that and let be as in Proposition 7.10

if then return

if , then

else ,

return

Proof. Let us first show that

is finite. Given

, we either have

and

, or

for some

. In the latter case, our assumption

implies the existence of a maximal

with

. We claim that there also exists a maximal

with

. Our claim clearly implies that

is finite, since we are left with a finite number of possible

for each

Assume for contradiction that the claim does not hold. Then

, whence there exists a

with

. For each

, there exists an

with

. For some

and

, it follows that

. Since

, we have

and

acts on

. In particular, there exist minimal

and maximal

with

. For

, we conclude that there cannot exist an

with

Having shown that

is finite, we either hit some

with

, in which case the correctness follows from Proposition 7.10, or all computations take place in the narrow submodule

, and the correctness and termination are proved in a similar way as for Theorem 6.4.□

8.Complexity of rational function arithmetic

8.1.Complexity of arithmetic in

Consider a standard DD-operator algebra

with

, as in section 6.1. We will use the algebraic complexity model: running times are measured in terms of the required number of field operations in

and space complexity in terms of the required number of coefficients in

We use a dense representation for polynomials in

. A polynomial

of total degree

thus requires space

. From now on, we assume that the dimension

is a fixed constant, whereas

may become large. Given

of total degrees

and

, it is well known for that their product can be computed in quasi-linear time

; see [53] for a particularly efficient algorithm. Partial derivatives and

-differences can clearly be computed in linear time

. It is also well-known [11] that partial shifts reduce to coefficientwise univariate multiplications, so they can be computed in time

The deterministic computation of the greatest common divisor of two polynomials

is more expensive. The best currently known algorithms first select a principal variable, say

, and reinterpret

and

as univariate polynomials in

. They next recursively factor out the content and then compute the univariate gcd using a division-free algorithm for the computation of subresultants such as Berkowitz' algorithm [8]; see also [69]. This leads to the complexity bound

with

. When using randomized algorithms of Las Vegas type, one may take

instead; since

, this means that gcds can be computed in quasi-linear time (more precisely: after a generic linear change of variables, both

and

can be assumed to be monic in

, and using evaluation/interpolation at

points, one reduces to fast univariate gcds [18]).

Rational functions in

are represented as fractions

with

, and

monic with respect to a suitable term ordering. We write

for the degree of such a rational function. Given

and

, we notice that

. One may compute

in non simplified form in quasi-linear time

, whereas the simplification of the resulting fraction requires

further operations. Of course, it is possible to delay the simplification of fractions in many cases.

When computing with vectors of matrices with entries in

, we always assume that denominators have been factored out. In other words, a matrix

is represented as a matrix

divided by a monic polynomial

such that the gcd of

and all entries of

equals one (here

is understood to be monic with respect to a suitable monomial ordering on

). We also define

. Simplifying a general fraction

can be done in time

. Ignoring the cost of the final simplification, two matrices

can be added in time

and multiplied in time

. Here

denotes the “exponent of matrix multiplication”, i.e. two matrices in

can be multiplied using

operations in

When allowing for probabilistic algorithms (of Las Vegas type), it is interesting to consider the alternative straight line program (SLP) representation for polynomials and rational functions [20]. In this model, the length

of the SLP that evaluates

is a natural mesure for the complexity of

. Using sparse interpolation [7, 58, 57, 51, 5], the actual coefficients of a polynomial

of degree bounded by

can be retrieved from its SLP representation in expected time

. Similarly, the dense representation of a rational function

can be retrieved in expected time

; see [59].

Algorithms that rely on SLP representations and sparse interpolation are particularly interesting when intermediate expressions swell in size, but the end-result is small in comparison. It turns out that the sizes of SLPs for creative telescoping remain reasonably small, which makes this approach relevant in our context. For this reason, we will also analyze the complexity of our algorithms in this model.

8.2.Complexity of arithmetic in

Now consider a standard DD-operator algebra extension

as in section 6.3. Polynomials in

and rational fractions in

can be represented in a similar way as in the previous subsection, but for our complexity analysis, it is useful to distinguish between the degree in

and the total degree in

. As in the case of vectors and matrices, we also assume that the denominators of polynomials in

have been factored out.

With the above conventions, the monic gcd

of two polynomials

can be computed in time

and we have

Similarly, given two polynomials

with monic

, the Euclidean division of

yields a relation

8.3.Complexity of arithmetic in

The computation of the local reductions

in sections 4 and 5 requires us to compute with roots

of the leading coefficient

of our differential or difference operator

. If the factorization of

is known, then this comes down to computations over the extension field

, where

is the monic irreducible factor of

with

Let us first examine the cost of arithmetic in such an extension field

and denote

. Elements

are represented by their pre-images

with

and

. It is natural to define the total degree of

. Given

, it follows from (8.6) that

where we recall from (8.2) that

. Similarly, the inverse of a non-zero

satisfies

thanks to (8.3). Since additions and subtraction do not require any Euclidean divisions, we also have

The computations of

and

can respectively be performed in time

and

, whereas

can be computed in time

. One may take

in the first two complexity bounds if we only require unsimplified results.

For simplicity, we will assume that the irreducible factorization of

is known in what follows. In order to be complete, we still have to discuss the cost of computing such a factorization. Since this task can be rather expensive, it is actually better to avoid it by relying on the strategy of “dynamic evaluation” [31].

Starting with the square-free part

, the idea is to directly work in the extension

as if it were a field. Whenever we hit a non-trivial zero-divisor

during the computation of an inverse in

, we simply abort all computations and resume with

and

in the role of

. This can happen at most

times, so all complexity bounds that will be proved in the sequel have to be multiplied by

in order to take into account the cost of the factorization of

. In fact, the

overhead can be reduced at the expense of technical complications; since we merely want to establish polynomial complexity bounds, we can spare ourselves this effort.

Computations with algebraic functions in

are often more efficient in the SLP model. In that case, at every evaluation point

, we work in the algebraic extension

instead of

. Here we used angular brackets for the evaluation of

in order to avoid confusion with evaluation in

. Now SLPs over this algebraic extension

can be reinterpreted as SLPs over

: additions/subtractions in

translate into

additions/subtractions in

, whereas each multiplication or division in

give rise to

ring operations in

. When implementing dynamic evaluation, one also must be careful to “take the same branch of the computation tree” at each evaluation point: whenever we hit a non-trivial zero-divisor

, we store the corresponding SLP for

, and first relaunch the computation with

in the role of

(for all evaluation points), and then once more with

in the role of

(again for all evaluation points).

8.4.Partial fraction decomposition

Consider a fraction

with

and

. Assume that

is monic and that its factorization over

is known, say

. Let us analyze the cost of computing the partial fraction decomposition of

The partial fraction decomposition can actually be computed over

as well as over its algebraic closure. When working over

, the partial fraction decomposition can be written at choice in one of the following two forms

where

and

. It is well known [65] that partial fraction decompositions of these kinds can be computed in quasi-linear time

in terms of operations over

In this paper, we have chosen to systematically work over the algebraic closure

. Now given a root

and any automorphism

over

, we notice that

. It therefore suffices to compute

for only one root

in each conjugacy class. We define the partial fraction decomposition of over to consist of

and the collection of coefficients

, where

runs over

and one element in each conjugacy class among the set of roots of

, and where

. The “tangling” and “untangling” morphisms from [52, Section 4] allow for conversions between partial fraction decompositions over

and

in quasi-linear time. It will be convenient to write

for the maximal degree in

of a coefficient in the partial fraction decomposition of

over

We have the following deterministic (but rather crude) bounds for the cost of partial fraction decomposition and the resulting degree swell in

Proof. From (8.5) and the fact that we may compute

as the quotient of the Euclidean division of

, we first notice that

for

. Let

be the monic minimal polynomial of

. We first notice that

for all

, whence

for all

and similarly for

. Let

and

be the power series expansions of

and

with

. It follows that

and

are bounded by

for all

, and

. Now consider the power series

. By what precedes, we have

and we notice that

for

. Since

and

, we conclude that

The quotient

of the Euclidean division of

and

can be computed in time

. The naive expansion of the series

and

until order

takes time

and the division of the results with precision

can be performed in time

. Doing this for one

in each conjugacy class yields the announced complexity bound.□

Partial fraction decomposition becomes much faster when carried out directly over

. In the SLP model, we have the following complexity bound:

Proof. Direct consequence of [65] and the quasi-linear time algorithm for the “untangling” morphism from [52, Sections 4.3 and 4.4].□

Since the local reductions really operate on partial fraction decompositions, it is also interesting to study the inverse operation that recovers a rational function from its partial fraction decomposition over

Proof. Let

be one the poles in the partial fraction decomposition of

and let

be its monic minimal polynomial. Denote

, and

. Consider the polar part

in the partial fraction decomposition of

with

and

. This components gives rise to a corresponding component

in the partial fraction decomposition of

over

. Now

, whence

and

. The traces of

can be computed using Newton–Girard's formulas and

for all

. Since

is linear, it follows that

Summing over all poles, the bound for

follows. The time complexity bound is proved in a similar way.□

Proof. We first convert the partial fraction decomposition over

back into a partial fraction decomposition over

; using the tangling algorithm from [52, Section 4.5], this can be done in quasi-linear time. We next recover the numerator and denominator using traditional rational function arithmetic, again in quasi-linear time.□

9.Complexity of creative telescoping

9.1.Complexity of creative telescoping

Let

be a DD-algebra as in the previous section. For all our complexity bounds, we recall that

is assumed to be constant. Consider a telescopable D-finite ideal

with associated reduction

. Let

with entries

be a basis for the finite dimensional vector space

. In Algorithm 6.1, it often occurs that around

iterations are necessary before the main loop terminates. In that case, it may be profitable to use the alternative way

to compute the

and to further optimize the reduction process by precomputing

for all basis elements

. For each

, this leads to a matrix

such that

Given

, we then get

. We call

the reduction matrices for

and

. With these optimizations, Algorithm 6.1 becomes:

Algorithm 9.1

Input: reduction matrices for and with

Output: D-finite generators of a D-finite ideal included in

, ,

while do

Let be minimal for and set

if is not reducible with respect to then

if then

else decompose and let

if , then

else ,

return

Remark 9.1. Finding sharp bounds for the complexity of the Gröbner basis computation in the final step of Algorithm 6.1 is a problem that is somewhat independent from the cost of creative telescoping itself. For this reason, we removed this step here. The returned set

still admits the property that the set of monomials that are not divisible by a leading monomial of an element of

is finite; for this reason we call

a D-finite set of generators.

Proof. Since each

is obtained as the product of at most

matrices of total degree

, we have

for all

encountered in the algorithm. The number of iterations of the main loop is always bounded by

. It follows that the computation of the

can be done in time

. The linear algebra of testing the existence of a relation

(and computing such a relation if it exists) reduces to the computation of the kernel of a matrix in

of total degree

and rank

, where

. This can be done in time

through evaluation-interpolation of the numerator at

points, using the fact that the determinants of all minors of this numerator have total degree

. Multiplying with the maximal number of iterations

, the main result follows.□

Proof. The evaluation of the input SLPs at a point

takes

operations and yields scalar matrices

and a scalar vector

. Applying the algorithm to

instead of

produces

instead of

. This time, each iteration amounts to one matrix-vector product of cost

and adding one row to the matrix of the

; we incrementally put this matrix in echelon normal form in order to determine the relations

, again with cost

. Since the algorithm uses at most

iterations, the complexity bound follows.□

Remark 9.4. The cost

of the linear algebra in Theorem 9.3 can be further lowered if the matrices

are sparse or in the case when we are using a lexicographic admissible ordering. For details, see [38] and references therein.

It may happen that some of the basis elements are actually superfluous for the computation, in which case it is possible to replace the reduction matrices

by smaller “quasi-reduction” matrices. More precisely, assume that the vector space

spanned by the basis elements

is stable under the mappings

for

. Then each

is a block matrix

where

. If

, then we call

quasi-reduction matrices for

and

. Given

with

, it is not hard to check that Algorithm 9.1 still applies to

instead of

, and that the above theorems generalize to this case.

9.2.Deterministic complexity of differential reduction

In this subsection, we assume that

and that our D-finite ideal

is explicitly telescopable for some cyclic basis

. We write

for the associated reduction. There also exists an operator

of order

than annihilates the first entry

. We may write

, where

and

is the monic annihilator of the zero-set of

. We assume that the irreducible factors of

are all known.

In order to apply Theorem 9.2, we need to bound the dimension

of the space

of reduced functions and to analyze the cost of computing the reduction matrices

. Let us start with the dimension bound.

Proof. We may take

in the proof of Proposition 7.5 (see Example 7.6). It follows that

. When applying Theorem 4.5 to

, we thus have

, whence the complexity bound, since

.□

Given a vector

with entries

, let us now analyze the cost to compute the reduction

. By construction, we recall that

, where

and

denotes the reduction with respect to the adjoint operator

. It thus suffices to analyze the cost of the reduction

Proof. Consider one of the poles

with

, let

be its monic minimal polynomial, and

. Given

, written as a Laurent polynomial, let

denote the maximal degree in

of one of its coefficients. Let us show by induction on

that

Let us now bound the quantity

. Rewriting

in terms of

and

amounts to a Taylor shift from which it is also possible to read off the indicial polynomial

. This yields

Let us now turn our attention to a general

. Proposition 8.1 implies that

. From what precedes, we have

at every pole

with

, whence

Proof. Denote

. For each

, let

be minimal such that

and

for

. Since

admits at most

roots and

, we notice that

. Consider the

-subvector space

and the corresponding

-subvector space

. The way we defined

ensures that we actually have

, whence

and

for

. This means that

is stable under the mappings

. We construct

as the quasi-reduction matrices with respect to a basis

The basis elements

can all be taken to be of the form

, where

is the minimal polynomial of some

and

. For each basis element

, it then follows that

. Writing

for the column vector with entries

, we also have

and

for

. Similarly, we have

, and

. Applying Lemma 9.6 to

, it follows that

. Regrouping the coefficients of the partial fraction decomposition of

by conjugate roots, we directly obtain the expression of

in terms of the basis elements

. Using Proposition 8.3, it follows that

We need to compute

for

and

, which means that we need

applications of Lemma 9.6, in view of Lemma 9.5. The complexity bound then follows from Lemma 9.6 and Proposition 8.3.□

Proof. In view of the discussion at the end of section 8.3, this is a direct consequence of Lemmas 9.5 and 9.7, together with the straightforward generalization of Theorem 9.2 to the case of quasi-reduction matrices.□

Remark 9.9. We did not put a lot of effort in lowering the exponents in the complexity bounds in Lemmas 9.6 and 9.7. One idea that should allow for significant improvements is to avoid conversions between the default representation of rational functions and partial fraction decompositions: it should be possible to carry out most computations directly for the partial fraction representation; see [36] for a similar line of thought in a different context.

9.3.SLP complexity of differential reduction

With the same assumptions as in the previous subsection, let us now carry out the complexity analysis in the SLP model. Before stating the counterpart of Lemma 9.6, let us first show that the local reductions from section 4.2 can essentially be computed in quasi-linear time.

Proof. The operator

can be rewritten as an operator in

using a Taylor shift; this can be done in time

. Rewritten in this way, the operator

naturally operates on Laurent polynomials in

. We may also decompose

with

and, for all

This decomposition can be done in linear time since both

and

are of the form

for suitable sets

. Now consider the linear operators

that act on monomials by

The second equation is “recursive”, which allows us to compute its solution in quasi-linear time in terms of the size of the equation using relaxed evaluation [48, 47]. This directly implies the first complexity bound. The bound at infinity is proved in a similar way.□

Proof. Denote

and let

be an evaluation point. By Proposition 8.2, we can compute the rational fraction decomposition of

over

in quasi-linear time

. For each

, let

be a root of

. Lemma 9.10 allows us to compute

in time

. Using the fact that

it follows that all

for

can be computed in time

. Using Proposition 8.4, these values allow us to compute

in quasi-linear time

. Again by Lemma 9.10, we obtain

using

more operations in

.□

Proof. Similarly as in the proof of Lemma 9.7, we apply Lemma 9.11 to

for

and

. We again have

, whence each individual

can be treated in time

. Since there are

values of

to consider, the conclusion follows.□

Proof. This is a direct consequence of Lemmas 9.5 and 9.12, together with the straightforward generalization of Theorem 9.3 to the case of quasi-reduction matrices.□

Remark 9.14. The bound

in the proof of Lemma 9.12 is actually quite pessimistic, since it corresponds to highly unlucky values of the roots of the indicial polynomials. In practice, the better bound

is likely to hold. In that case, the time complexity per evaluation point in Theorem 9.13 drops to

Bibliography

[1]: S. A. Abramov. Eg-eliminations. Journal of Difference Equations and Applications, 5(4–5):393–433, 1999.
[2]: S. A. Abramov and M. Petrovšek. Rational normal forms and minimal decompositions of hypergeometric terms. JSC, 33(5):521–543, 2002.
[3]: K. Adjamagbo. Sur l'effectivité du lemme du vecteur cyclique. C. R. Acad. Sci. Paris Sér. I, 306(13):543–546, 1988.
[4]: G. Almkvist and D. Zeilberger. The method of differentiating under the integral sign. Journal of Symbolic Computation, 10(6):571–591, 1990.
[5]: A. Arnold, M. Giesbrecht, and D. S. Roche. Faster sparse multivariate polynomial interpolation of straight-line programs. Technical Report http://arxiv.org/abs/1412.4088, Arxiv, 2015. To appear in JSC.
[6]: M. Aschenbrenner, L. van den Dries, and J. van der Hoeven. Asymptotic Differential Algebra and Model Theory of Transseries. Number 195 in Annals of Mathematics studies. Princeton University Press, 2017.
[7]: M. Ben-Or and P. Tiwari. A deterministic algorithm for sparse multivariate polynomial interpolation. In STOC '88: Proceedings of the twentieth annual ACM symposium on Theory of computing, pages 301–309. New York, NY, USA, 1988. ACM Press.
[8]: S. J. Berkowitz. On computing the determinant in small parallel time using a small number of processors. Inform. Process. Lett., 18:147–140, 1984.
[9]: I. N. Bernshtein. Modules over a ring of differential operators, study of the fundamental solutions of equations with constant coefficients. Functional Anal. Appl., 5(2):89–101, 1971.
[10]: I. N. Bernshtein. The analytic continuation of generalized functions with respect to a parameter. Functional Anal. Appl., 6(4):273–285, 1972.
[11]: D. Bini and V. Y. Pan. Polynomial and matrix computations. Vol. 1. Birkhäuser Boston Inc., Boston, MA, 1994. Fundamental algorithms.
[12]: J. E. Björk. Rings of differential operators. North-Holland, New York, 1979.
[13]: A. Bostan, F. Chen, S. Chyzak, and Z. Li. Complexity of creative telescoping for bivariate rational functions. In Proc. ISSAC '12, pages 203–210. New York, NY, USA, 2010. ACM.
[14]: A. Bostan, S. Chen, F. Chyzak, Z. Li, and G. Xin. Hermite reduction and creative telescoping for hyperexponential functions. In Proc. ISSAC'13, pages 77–84. ACM, 2013.
[15]: A. Bostan, F. Chyzak, and É. de Panafieu. Complexity estimates for two uncoupling algorithms. In Proc. ISSAC 2013, ISSAC '13, pages 85–92. New York, NY, USA, 2013. ACM.
[16]: A. Bostan, F. Chyzak, P. Lairez, and B. Salvy. Generalized Hermite reduction, creative telescoping and definite integration of differentially finite functions. Technical Report, ArXiv, 2018. http://arxiv.org/abs/1805.03445.
[17]: F. Boulier. Étude et implantation de quelques algorithmes en algèbre différentielle. PhD thesis, University of Lille I, 1994.
[18]: R. P. Brent, F. G. Gustavson, and D. Y. Y. Yun. Fast solution of Toeplitz systems of equations and computation of Padé approximants. J. Algorithms, 1(3):259–295, 1980.
[19]: B. Buchberger. Ein Algorithmus zum auffinden der Basiselemente des Restklassenringes nach einem null-dimensionalen Polynomideal. PhD thesis, University of Innsbruck, 1965.
[20]: Bürgisser, P. and Clausen, M. and Shokrollahi, M. A. Algebraic complexity theory, volume 315 of Grundlehren der Mathematischen Wissenschaften. Springer-Verlag, 1997.
[21]: S. Chen, M. van Hoeij, M. Kauers, and C. Koutschan. Reduction-based creative telescoping for Fuchsian D-finite functions. Technical Report, ArXiv, 2016. http://arxiv.org/abs/1611.07421.
[22]: S. Chen, H. Huang, M. Kauers, and Z. Li. A modified Abramov-Petkovs^ˇek reduction and creative telescoping for hypergeometric terms. In Proc. ISSAC'15, pages 117–124. ACM, 2015.
[23]: S. Chen, M. Kauers, and C. Koutschan. Reduction-based creative telescoping for algebraic functions. In Proc. ISSAC '16, pages 175–182. New York, NY, USA, 2016. ACM.
[24]: S. Chen, M. Kauers, and M. Singer. Telescopers for rational and algebraic functions via residues. In Proc. ISSAC'12, pages 130–137. ACM, 2012.
[25]: R. C. Churchill and J. J. Kovacic. Differential Algebra and Related Topics, chapter Cyclic vectors, pages 191–218. World Scientific, 2002.
[26]: F. Chyzak. Fonctions holonomes en calcul formel. PhD thesis, École polytechnique, France, 1998.
[27]: F. Chyzak. An extension of Zeilberger's fast algorithm to general holonomic functions. Discrete Math., 217(1–3):115–134, 2000.
[28]: F. Chyzak. The ABC of Creative Telescoping — Algorithms, Bounds, Complexity. Habilitation, École polytechnique, 2014.
[29]: F. Chyzak and B. Salvy. Non–commutative elimination in Ore algebras proves multivariate identities. JSC, 26(2):187–227, 1998.
[30]: J. H. Davenport. The Risch differential equation problem. SIAM J. Comput., 15(4):903–918, 1986.
[31]: J. Della Dora, C. Dicrescenzo, and D. Duval. A new method for computing in algebraic number fields. In G. Goos and J. Hartmanis, editors, Eurocal'85 (2), volume 174 of Lect. Notes in Comp. Science, pages 321–326. Springer, 1985.
[32]: J. Denef and L. Lipshitz. Power series solutions of algebraic differential equations. Math. Ann., 267:213–238, 1984.
[33]: J. Denef and L. Lipshitz. Decision problems for differential equations. The Journ. of Symb. Logic, 54(3):941–950, 1989.
[34]: L. Dumont. Efficient algorithms for the symbolic computation of some contour integrals depending on one parameter. PhD thesis, École Polytechnique, 2016.
[35]: M. C. Fasenmyer. Some generalized hypergeometric polynomials. PhD thesis, Univ. of Michigan, 1945.
[36]: R. J. Fateman. Rational function computing with poles and residues. http://www.cs.berkeley.edu/~fateman/papers/poles.pdf, 2010.
[37]: J. C. Faugère, P. Gianni, D. Lazard, and T. Mora. Efficient computation of zero-dimensional Gröbner bases by change of ordering. JSC, 16(4):329–344, 1993.
[38]: J.-C. Faugère, P. Gaudry, L. Huot, and G. Renault. Sub-cubic Change of Ordering for Gröbner Basis: A Probabilistic Approach. In Proc. ISSAC 2014, pages 170–177. Kobe, Japon, jul 2014. ACM.
[39]: L. Fuchs. Die Periodicitaätsmoduln der hyperelliptischen Integrale als Functionen eines Parameters aufgefasst. J. Reine Angew. Math., 71:91–127, 1870.
[40]: A. Galligo. Some algorithmic questions on ideals of differential operators. In Proc. EUROCAL '85, volume 204 of Lecture Notes in Computer Science, pages 413–421. Springer-Verlag, 1985.
[41]: K. Geddes, H. Le, and Z. Li. Differential rational normal forms and a reduction algorithm for hyperexponential functions. In Proc. ISSAC'04, pages 183–190. ACM, 2004.
[42]: R. W. Gosper, Jr. Decision procedure for indefinite hypergeometric summation. Proc. Nat. Acad. Sci. U.S.A., 75(1):40–42, 1978.
[43]: G.-M. Greuel, V. Levandovskyy, O. Motsak, and H. Schönemann. Plural. a singular 3-1 subsystem for computations with non-commutative polynomial algebras. http://www.singular.uni-kl.de, 2010.
[44]: C. Hermite. Sur l'intégration des fractions rationnelles. Ann. Sci. École Norm. Sup. Série 2, 1:215–218, 1972.
[45]: J. van der Hoeven. Differential and mixed differential-difference equations from the effective viewpoint. Technical Report LIX/RR/96/11, LIX, École polytechnique, France, 1996.
[46]: J. van der Hoeven. A new zero-test for formal power series. In Teo Mora, editor, Proc. ISSAC '02, pages 117–122. Lille, France, July 2002.
[47]: J. van der Hoeven. Relax, but don't be too lazy. JSC, 34:479–542, 2002.
[48]: J. van der Hoeven. New algorithms for relaxed multiplication. JSC, 42(8):792–802, 2007.
[49]: J. van der Hoeven. Computing with D-algebraic power series. Technical Report, HAL, 2014. http://hal.archives-ouvertes.fr/hal-00979367, accepted for publication in AAECC.
[50]: J. van der Hoeven. Constructing reductions for creative telescoping. Technical Report, HAL, 2017. http://hal.archives-ouvertes.fr/hal-01435877.
[51]: J. van der Hoeven and G. Lecerf. Sparse polynomial interpolation in practice. ACM Commun. Comput. Algebra, 48(3/4):187–191, 2015.
[52]: J. van der Hoeven and G. Lecerf. Composition modulo powers of polynomials. In Proc. ISSAC '17, pages 445–452. New York, NY, USA, 2017. ACM.
[53]: J. van der Hoeven and É. Schost. Multi-point evaluation in higher dimensions. AAECC, 24(1):37–52, 2013.
[54]: H. Huang. New bounds for hypergeometric creative telescoping. In Proc. ISSAC'16, pages 279–286. ACM, 2016.
[55]: E. L. Ince. Ordinary differential equations. Longmans, Green and Co., 1926. Reprinted by Dover in 1944 and 1956.
[56]: M. Janet. Sur les systèmes d'équations aux dérivées partielles. PhD thesis, Faculté des sciences de Paris, 1920. Thèses françaises de l'entre-deux-guerres. Gauthiers-Villars.
[57]: M. Javadi and M. Monagan. Parallel sparse polynomial interpolation over finite fields. In Proceedings of PASCO 2010, pages 160–168. ACM Press, 2010.
[58]: E. Kaltofen, Y. N. Lakshman, and J.-M. Wiley. Modular rational sparse multivariate polynomial interpolation. In ISSAC '90: Proceedings of the international symposium on Symbolic and algebraic computation, pages 135–139. New York, NY, USA, 1990. ACM Press.
[59]: E. Kaltofen and B. M. Trager. Computing with polynomials given by black boxes for their evaluations: greatest common divisors, factorization, separation of numerators and denominators. JSC, 9(3):301–320, 1990.
[60]: A. Kandri-Rody and V. Weispfenning. Non-commutative Gröbner bases in algebras of solvable type. JSC, 9:1–26, 1990.
[61]: E. R. Kolchin. Differential algebra and algebraic groups. Academic Press, New York, 1973.
[62]: T. H. Koornwinder. On Zeilberger's algorithm and its q-analogue. In Proceedings of the Seventh Spanish Symposium on Orthogonal Polynomials and Applications (VII SPOA), volume 48, pages 91–111. Granada (1991), 1993.
[63]: C. Koutschan. Advanced Applications of the Holonomic Systems Approach. PhD thesis, RISC-Linz, 2009.
[64]: H. Kredel. Solvable polynomial rings. PhD thesis, Univ. Passau, 1993. Published by Shaker.
[65]: H. T. Kung and D. M. Tong. Fast algorithms for partial fraction decomposition. SIAM J. Comput., 6:582–593, 1977.
[66]: J.-L. Lagrange. Miscellanea taurinensia iii. p. 179.
[67]: P. Lairez. Periods of rational integrals: algorithms and applications. PhD thesis, École polytechnique, Nov 2014.
[68]: P. Lairez. Computing periods of rational integrals. Math. Comp., 85(300):1719–1752, 2016.
[69]: G. Lecerf. On the complexity of the Lickteig-Roy subresultant algorithm. Technical Report, CNRS & École polyechnique, 2017. https://hal.archives-ouvertes.fr/hal-01450869.
[70]: V. Levandovskyy. Non-commutative Computer Algebra for polynomial algebras: Gröbner bases, applications and implementation. PhD thesis, Univ. Kaiserslautern, 2005.
[71]: T. Mora. Groebner bases in non-commutative algebras. In Proc. ISSAC'88, number 358 in LNCS, pages 150–161. 1989.
[72]: O. Ore. Theorie der linearen Differentialgleichungen. J. Reine und Angew. Math., 167:221–234, 1932.
[73]: O. Ore. Theory of non-commutative polynomials. Ann. of Math., 34(3):480–508, 1933.
[74]: M. Ostrogradsky. De l'integration des fractions rationelles. Bull. de la Classe Physico- Mathématique de l'Académie Imperiale des Sciences de Saint-Petersburg, IV:147–168, 1845.
[75]: A. Péladan-Germa. Tests effectifs de nullité dans des extensions d'anneaux différentiels. PhD thesis, Gage, École Polytechnique, Palaiseau, France, 1997.
[76]: É. Picard. Sur les intégrales doubles de fonctions rationnelles dont tous les résidus sont nuls. Bull. Sci. Math. (Série 2), 26:143–152, 1902.
[77]: A. van der Poorten. A proof that Euler missed: Apéry's proof of the irrationality of . Math. Intelligencer, 1(4):195–203, 1979.
[78]: C. G. Raab. Definite Integration in Differential Fields. PhD thesis, Johannes Kepler Universität Linz, 2012.
[79]: Ch. Riquier. Les systèmes d'équations aux dérivées partielles. Gauthier-Villars, 1910.
[80]: J. F. Ritt. Differential algebra. Amer. Math. Soc., New York, 1950.
[81]: A. Rosenfeld. Specializations in differential algebra. Trans. Amer. Math. Soc., 90:394–407, 1959.
[82]: C. Schneider. Simplifying Multiple Sums in Difference Fields. In C. Schneider and J. Bluemlein, editors, Computer Algebra in Quantum Field Theory: Integration, Summation and Special Functions, Texts and Monographs in Symbolic Computation, pages 325–360. Springer, 2013.
[83]: J. Shackell. A differential-equations approach to functional equivalence. In Proc. ISSAC '89, pages 7–10. Portland, Oregon, A.C.M., New York, 1989. ACM Press.
[84]: N. Takayama. Gröbner basis and the problem of contiguous relations. Japan J. Appl. Math., 6:147–160, 1989.
[85]: N. Takayama. Kan: a system for computation in algebraic analysis. http://www.math.kobe-u.ac.jp/KAN/, 1991.
[86]: N. Takayama. An approach to the zero recognition problem by Buchberger algorithm. JSC, 14:265–282, 1992.
[87]: J. M. Thomas. Differential Systems, volume XXI. American Mathematical Society Colloquium Publications, New York, 1937.
[88]: B. M. Trager. Integration of Algebraic Functions. PhD thesis, MIT, 1984.
[89]: P. Verbaeten. Rekursiebetrekkingen voor lineaire hypergeometrische funkties. PhD thesis, K. U. Leuven, 1976.
[90]: D. Zeilberger. A holonomic systems approach to special functions identities. Journal of Comp. and Appl. Math., 32:321–368, 1990.
[91]: D. Zeilberger. The method of creative telescoping. JSC, 11(3):195–204, 1991.

			(θ_i=∂_i,θ_j=∂_j)
			(θ_i=σ_i,θ_j=∂_j)
			(θ_i=σ_i,θ_j=σ_j)

			(2.4)
			(2.5)

			(α∈𝕂^¯)
			(α=∞)

			(α∈𝕂^¯)
			(α=∞)

1.Introduction

1.1.Creative telescoping

1.2.From Hermite reduction to reduction-based telescoping

1.3.On the choice of our setup

1.4.Structure of the paper and outline of the main results

1.5.Pending issues

2.Abstract reductions

2.1.Definitions and basic properties

2.2.Local reductions

2.3.Gluing local reductions together

2.4.Composition of reductions

2.5.Normalization of confined reductions

3.Lagrange's identity and generalizations

3.1.The original differential case

3.2.The difference case

3.3.Matrix versions

3.4.Twisting

3.5.Another twisted variant

4.Differential reductions

4.1.Partial fraction decomposition and Hermite reduction

4.2.Local reduction with respect to differential operators

4.3.Global reduction with respect to differential operators

4.4.Back to reductions with respect to

4.5.First order systems

4.6.Reductions on D-modules with respect to general operators

5.Difference reductions

5.1.The field as a difference field

5.2.Difference modules

5.3.Local and semi-local reduction with respect to shift operators

5.4.Aligned downward reduction with respect to shift operators

5.5.General downward reduction with respect to shift operators

5.6.Back to reductions with respect to

5.7.The case of -difference operators

5.8.First order systems

6.Creative telescoping

6.1.D-finite ideals of DD-operator algebras

6.2.Matrix representations of operators

6.3.Narrowness

6.4.Creative telescoping

7.D-finiteness tests for telescoping ideals

7.1.Location of singularities

7.2.Orders at singularities

7.3.Explicit telescopability

7.4.Testing D-finiteness

8.Complexity of rational function arithmetic

8.1.Complexity of arithmetic in

8.2.Complexity of arithmetic in

8.3.Complexity of arithmetic in

8.4.Partial fraction decomposition

9.Complexity of creative telescoping

9.1.Complexity of creative telescoping

9.2.Deterministic complexity of differential reduction

9.3.SLP complexity of differential reduction

Bibliography