Hermite Reduction and Creative Telescoping for Hyperexponential Functions
Alin Bostan, Shaoshi Chen, Frédéric Chyzak, Ziming Li, Guoce Xin
aa r X i v : . [ c s . S C ] J a n Hermite Reduction and Creative Telescopingfor Hyperexponential Functions ∗ Alin Bostan , Shaoshi Chen , Frédéric Chyzak , Ziming Li , Guoce Xin INRIA, Palaiseau, 91120, (France) Department of Mathematics, NCSU, Raleigh, 27695-8025, (USA) KLMM, AMSS, Chinese Academy of Sciences, Beijing 100190, (China) Department of Mathematics, Capital Normal University, Beijing 100048, (China) [email protected], { alin.bostan, frederic.chyzak } @[email protected], [email protected] ABSTRACT
We present a reduction algorithm that simultaneously ex-tends Hermite’s reduction for rational functions and theHermite-like reduction for hyperexponential functions. Ityields a unique additive decomposition and allows to de-cide hyperexponential integrability. Based on this reductionalgorithm, we design a new method to compute minimaltelescopers for bivariate hyperexponential functions. One ofits main features is that it can avoid the costly computa-tion of certificates. Its implementation outperforms Maple’sfunction
DEtools[Zeilberger] . Moreover, we derive an orderbound on minimal telescopers, which is more general andtighter than the known one.
Categories and Subject Descriptors
I.1.2 [
Computing Methodologies ]: Symbolic and Alge-braic Manipulation—
Algebraic Algorithms
General Terms
Algorithms, Theory
Keywords
Hermite Reduction, Hyperexponential function, Telescoper
1. INTRODUCTION
Given a univariate rational function r , Hermite reduc-tion in [12, 13, 6] finds rational functions r and r s.t. ∗ S.C. was supported by the National Science Foundation(NSF) grant CCF-1017217, A.B. and F.C. were supportedin part by the MSR-INRIA Joint Centre, Z.L. by two NSFCgrants (91118001 and 60821002/F02), and G.X. by NSFCgrant (11171231).
Permission to make digital or hard copies of all or part of this work forpersonal or classroom use is granted without fee provided that copies arenot made or distributed for profit or commercial advantage and that copiesbear this notice and the full citation on the first page. To copy otherwise, torepublish, to post on servers or to redistribute to lists, requires prior specificpermission and/or a fee.Copyright 20XX ACM X-XXXXX-XX-X/XX/XX ...$10.00. (i) r = r + r , (ii) r is rational integrable, (iii) r is aproper fraction with a squarefree denominator. The addi-tive decomposition is unique, and r is rational integrable ifand only if r = 0.A univariate function is hyperexponential if its logarith-mic derivative is rational. Exponential, radical and rationalfunctions are hyperexponential. Rational Hermite reductionhas been extended to hyperexponential functions by Daven-port in [9] and by Geddes, Le and Li in [10]. The formeraims at solving Risch’s equation; the latter is a differen-tial analogue of the reduction algorithm for hypergeometricterms in [2]. For a given hyperexponential function H , thereduction algorithms in [9, 10] compute two hyperexponen-tial functions H and H s.t. (i) H = H + H , (ii) H is hyperexponential integrable, (iii) H is minimal in somesense. However, H is not unique in general and it may benonzero even when H is hyperexponential integrable. In or-der to decide the integrability of H , one additionally needsto compute polynomial solutions of a first-order linear dif-ferential equation.The method of creative telescoping for hyperexponentialfunctions is developed by Almkvist and Zeilberger in [3]. Itis nowadays an important automatic tool for computing def-inite integrals. Recently, it has also played an important rolein the resolution of intriguing problems in enumerative com-binatorics [14, 15]. For a bivariate hyperexponential func-tion H ( x, y ), the problem of creative telescoping is to find anonzero operator L ( x, D x ) ∈ F ( x ) h D x i , the ring of linear dif-ferential operators over the rational-function field F ( x ), s.t. L ( x, D x )( H ) = D y ( G ) (1)for some hyperexponential function G , where D x = ∂/∂x and D y = ∂/∂y . The operator L above is called a telescoper for H , and G is the corresponding certificate . An algorithmfor solving (1) is given in [3], and is based on differentialGosper’s algorithm. An algorithm for rational-function tele-scoping is given in [5], and is based on Hermite reduction.The latter separates the computation for telescopers fromthat for certificates, and has a lower complexity than theformer for rational functions.In the present paper, we develop a reduction algorithmwhich, given a univariate hyperexponential function H , con-structs two hyperexponential functions H and H s.t. (i) H = H + H , (ii) H is hyperexponential integrable, andiii) H is either zero or not hyperexponential integrable. Weshow that H in the above additive decomposition is uniqueand can be obtained without computing polynomial solu-tions of any differential equation. Our algorithm is basedon the Hermite-like reduction in [10], a differential variantof the polynomial reduction in [2] and on the idea for re-ducing simple radicals in [18, Proposition 7]. The main newingredients are the uniqueness of H and an easy way tocompute H , which are crucial for many applications. Theseenable us to extend the reduction-based rational telescopingalgorithm in [5] to the hyperexponential case, and derive anorder bound on the telescopers. The bound is more generaland tighter than that given in [4].The rest of the paper is organized as follows. We reviewthe notion of hyperexponential functions and Hermite-likereduction in Sections 2 and 3, respectively. A new reduc-tion algorithm is developed for hyperexponential functionsin Section 4. After introducing kernel reduction in Sec-tion 5, we present a reduction-based telescoping algorithmfor bivariate hyperexponential functions, and derive an up-per bound on the order of minimal telescopers in Section 6.We briefly describe an implementation of the new telescop-ing algorithm, and present some experimental results in Sec-tion 7, which validate its practical relevance.As a matter of notation, we let E be a field of characteristiczero and E ( y ) be the field of rational functions in y over E .For a polynomial p ∈ E [ y ], we denote by deg( p ) and lc( p )the degree and leading coefficient of p , respectively. Let D y denote the usual derivation d/dy on E ( y ). Then ( E ( y ) , D y )is a differential field.
2. HYPEREXPONENTIAL FUNCTIONS
Hyperexponential functions share the common propertiesof rational functions, simple radicals, and exponential func-tions. Together with hypergeometric terms, they are fre-quently viewed as a special and important class of “closed-form” solutions of linear differential and difference equationswith polynomial coefficients.
Definition 1.
Let Φ be a differential field extension of E ( y ) .A nonzero element H ∈ Φ is said to be hyperexponential over E ( y ) if its logarithmic derivative D y ( H ) /H ∈ E ( y ) . The product of hyperexponential functions is also hyper-exponential. Two hyperexponential functions H , H aresaid to be similar if there exists r ∈ E ( y ) s.t. H = rH .The sum of similar hyperexponential functions is still hy-perexponential, provided that it is nonzero.For brevity, we use the notation exp( R fdy ) to indicate ahyperexponential function whose logarithmic derivative is f .For a rational function r ∈ E , we have r exp (cid:18)Z f dy (cid:19) = exp (cid:18)Z ( f + D y ( r ) /r ) dy (cid:19) . A univariate hyperexponential function H is said to be hy-perexponential integrable if it is the derivative of another hy-perexponential function. For brevity, we say “integrable” in-stead of “hyperexponential integrable” in the sequel.Assume that H = r exp (cid:0)R fdy (cid:1) is integrable. Then H is equal to D y ( G ) for some hyperexponential function G . Astraightforward calculation shows that G is similar to D y ( G ),and so is H . Set G = s exp (cid:0)R f (cid:1) for some s ∈ E ( y ).Then H = D y ( G ) if and only if r = D y ( s ) + f s. (2) Deciding the integrability of H amounts to finding a rationalsolution s s.t. the above equation holds.
3. HERMITE-LIKE REDUCTION
Reduction algorithms have been developed for computingadditive decompositions of rational functions [16, 12, 13], hy-pergeometric terms [1, 2], and hyperexponential functions [9,10]. Those algorithms can be viewed as generalizations ofGosper’s algorithm [11] and its differential analogue [3, § H , a reduction algorithmcomputes two hyperexponential functions H , H s.t. H = D y ( H ) + H . (3)It turns out that H, H and H are similar. So we maywrite H = r exp (cid:0)R fdy (cid:1) and H i = r i exp (cid:0)R fdy (cid:1) , where r, r i , f belong to E ( y ) and i = 1 ,
2. Then (3) translates into r = D y ( r ) + f r + r . A reduction algorithm for computing (3) amounts to choos-ing rational functions r, f and r so that r satisfies prop-erties similar to those obtained in Hermite reduction for ra-tional functions. There are at least two approaches to thisend. One is given in [9], and the other in [10]. We reviewthe latter, because the notion of differential-reduced rationalfunctions plays a key role in Lemma 6 in Section 4.Recall [10, §
2] that a rational function r = a/b ∈ E ( y ) issaid to be differential-reduced w.r.t. y ifgcd ( b, a − i D y ( b )) = 1 for all i ∈ Z .By Lemma 2 in [10], r is differential-reduced if and onlyif none of its residues is an integer. The differential ra-tional canonical form of a rational function f in E ( y ) is apair ( K, S ) in E ( y ) × E ( y ) s.t. (i) K is differential-reduced;(ii) the denominator of S is coprime with that of K ; and(iii) f is equal to K + D y ( S ) /S . Every rational function hasa unique canonical form in the sense that K is unique and S is unique up to a multiplicative constant in E [10, § K and S the kernel and shell of f , respectively. Theycan be constructed by the method described in [10, § H be a univariate hyperexponential function in theform exp( R fdy ) over E . Assume that K and S are the kerneland shell of f , respectively. Then H = S exp (cid:0)R K dy (cid:1) . Notethat K = 0 if and only if H is a rational function, which isequal to cS for some c ∈ E . Example 2.
Let H = p y + 1 / ( y − . The logarithmicderivative of H is D y HH = D y (1 / ( y − )1 / ( y − + yy + 1 , where y/ ( y +1) is differential-reduced. The kernel and shellof D y ( H ) /H are y/ ( y + 1) and / ( y − , respectively.So H = exp (cid:0)R y/ ( y + 1) dy (cid:1) / ( y − . For brevity, we make a notational convention.
Convention 3.
Let H denote a hyperexponential functionwhose logarithmic derivative has kernel K and shell S . As-sume that K is nonzero, that is, H is not a rational func-tion. Set T = exp (cid:0)R K dy (cid:1) . Moreover, write K = k /k ,where k , k are polynomials in E [ y ] with gcd( k , k ) = 1 . he algorithm ReduceCert in [10] computes a rationalfunction S s.t. S = D y ( S ) + S K + abk , (4)where a, b ∈ E [ y ] satisfy the following conditions: b is thesquarefree part of the denominator of S , and gcd( b, k )=1.Note that a is not necessarily coprime with bk . As the al-gorithm ReduceCert only reduces the shell S , it is referredto as the shell reduction . It follows from (4) that H = D y ( S T ) + abk T. (5)By Theorem 4 in [10], a/b belongs to E [ y ] if H is integrable. Example 4.
Let H be the same hyperexponential functionas in Example 2. Then D y ( H ) /H has kernel K = y/ ( y +1) and shell S = 1 / ( y − . The shell reduction yields S = D y ( S ) + S K + y ( y − k , where S = 1 / (1 − y ) and k = y +1 . Then H can be decom-posed into H = D y ( S T )+ yT / (( y − k ) , where T = p y +1 .By Theorem 4 in [10], H is not integrable. On the other hand, it is possible that a in (5) is nonzerobut H is integrable. Example 5.
Let H = y exp( y ) whose logarithmic derivativehas kernel and shell y , that is, H = y exp (cid:0)R dy (cid:1) . But H is integrable as it is equal to D y ( y exp( y ) − exp( y )) . The shell reduction cannot be directly used to decide hy-perexponential integrability. To amend this, the solutionproposed in [10, Algorithm
ReduceHyperexp ] was to findthe polynomial solutions of an auxiliary first-order lineardifferential equation. In the following section, we show howthis can be avoided and improved.
4. HERMITE REDUCTION FORHYPEREXPONENTIAL FUNCTIONS
After the shell reduction described in (5), the denomina-tors of shells have been reduced to squarefree polynomials.In the rational case, i.e., when the kernel K is zero, the poly-nomial a in (4) can be chosen s.t. deg( a ) < deg( b ), becauseall polynomials are rational integrable. But a hyperexpo-nential function with a polynomial shell is not necessarilyintegrable. For example, H = exp (cid:0) y (cid:1) . We present a differential variant of [2, Theorem 7] tobound the degree of a in (4). The variant leads not onlyto a canonical additive decomposition of hyperexponentialfunctions, but also a direct way to decide their integrability. With Convention 3, we define M K = { k D y ( p ) + k p | p ∈ E [ y ] } . It is an E -linear subspace in E [ y ]. We call M K the subspacefor polynomial reduction w.r.t. K . Moreover, define an E -linear map φ K from E [ y ] to M K that, for every p ∈ E [ y ],maps p to k D y ( p )+ k p . We call φ K the map for polynomialreduction w.r.t. K .Concerning the subspace M K and the map φ K , we have Lemma 6. (i) If k D y ( g ) + k g ∈ E [ y ] for some g ∈ E ( y ) , then g ∈ E [ y ] . (ii) The map φ K is bijective.Proof. Assume that g has a pole. Without loss of generality,we assume that the pole is y = 0 and has order m , becausethe following argument is also applicable over the algebraicclosure of E . Expanding g around the origin yields g = ry m + terms of higher orders in y, where r ∈ E \ { } . It follows from k D y ( g ) + k g ∈ E [ y ] that y = 0 is a pole of (cid:18) − mry m +1 + higher terms (cid:19) + K (cid:18) ry m + higher terms (cid:19) with order no more than that of K . This implies that y =0is a simple pole of K with residue m , which is incompatiblewith K being differential-reduced. The first assertion holds.The map φ K is surjective by its definition. If φ K ( p )=0 forsome nonzero polynomial p ∈ E [ y ], then K equals − D y ( p ) /p ,which is nonzero since K = 0. So K is not differential-reduced, a contradiction. The second assertion holds.An E -basis of M K is called an echelon basis if distinctelements in the basis have distinct degrees. Echelon basesalways exist and their degrees form a unique subset of N .Let B be an echelon basis of M K . Define N K = span E n x ℓ | ℓ ∈ N and ℓ = deg( f ) for all f ∈ B o . Then E [ y ] = M K ⊕ N K . We call N k the standard com-plement of M K . Using an echelon basis of M K , one canreduce a polynomial p to a unique polynomial ˜ p ∈ N K s.t. p − ˜ p ∈ M K .In order to find an echelon basis of M K , we set d = deg k , d = deg k , τ K = − lc( k ) / lc( k ), and B = { φ K ( y n ) | n ∈ N } .By Lemma 6 (ii), B is an E -basis of M K . Let p be a nonzeropolynomial in E [ y ]. We make the following case distinction. Case 1. d ≥ d . Then φ K ( p ) = lc( k ) lc( p ) y d +deg p + lower terms.So B is an echelon basis, in which deg φ K ( y n ) = d + n forall n ∈ N . Accordingly, N K is spanned by 1 , y, . . . , y d − . Case 2. d = d − τ K is not a positive integer. Then φ K ( p )= (deg( p ) lc( k )+ lc( k )) lc( p ) y d +deg p +lower terms . (6)Since τ K is not a positive integer, deg φ K ( y n ) = d + n .Thus, M K and N K have the same bases as in Case 1. Case 3. d < d −
1. If deg( p ) >
0, then φ K ( p ) = deg( p ) lc( k ) lc( p ) y d +deg( p ) − + lower terms . Otherwise, deg p = 0 and φ K ( p ) = k p . Therefore, B isagain an echelon basis, in whichdeg φ K (1) = d and deg φ K ( y n ) = d + n − n ≥ . Accordingly, N K has a basis 1 , . . . , y d − , y d +1 , . . . , y d − . Case 4. d = d − τ K is a positive integer. It followsfrom (6) that deg φ ( y n ) = d + n if n = τ K . Furthermore,for every polynomial p of degree τ K , deg( φ K ( p )) is of degreeless than d + τ K . So any echelon basis of M K does notcontain a polynomial of degree d + τ K . Set B ′ = { φ ( y n ) | n ∈ N , n = τ K } . educing φ ( y τ K ) by the polynomials in B ′ , we obtain apolynomial r of degree less than d . Note that r is nonzero,because B is an E -linearly independent set. Hence, B ′ ∪ { r } is an echelon basis of M K . Consequently, N K has an E -basis n , y, . . . , y deg( r ) − , y deg( r )+1 , . . . , y d − , y d + τ K o . Example 7.
Let K = − y / ( y + 1) , which is differential-reduced. Then τ K = 6 . According to Case 4, M K has anechelon basis (cid:8) y }∪{ ( n − y n +3 + ny n − | n ∈ N , n =6 (cid:9) . More-over, N K has a basis { , y , y } . One can reduce the degree and number of terms of a poly-nomial using a subspace of polynomial reduction.
Lemma 8.
With Convention 3, we further let d = deg k , d = deg k , and τ K = − lc( k ) / lc( k ) . Let M K be the sub-space for polynomial reduction, and N K its standard com-plement w.r.t. K . Finally, let p be a polynomial in E [ y ] .(i) If d ≥ d or d = d − and τ K / ∈ Z + , then thereexists q ∈ N K s.t. p ≡ q mod M K and deg q < d .(ii) If d < d − , then there exists q ∈ N K s.t. p ≡ q mod M K , deg q < d and the coefficient of y d in q is equal to zero.(iii) If d = d − and τ K ∈ Z + , then there exists r ∈ E [ y ] ofdegree less than d s.t. p ≡ sy d + τ K + r mod M K for some s ∈ E .Moreover, sy d + τ K + r belongs to N K , and r has atmost d − terms.Proof. The lemma is immediate from the E -bases of N K constructed in the above case distinction.The next corollary enables us to derive an order bound ontelescopers for hyperexponential functions. Corollary 9.
With the notation introduced in Lemma 8,there exists
P ⊂ { y n | n ∈ N } with |P| ≤ max( d , d − s.t.every polynomial in E [ y ] can be reduced modulo M K to an E -linear combination of the elements in P .Proof. By the above case distinction, the dimension of N K over E is at most max( d , d − With Convention 3, we further assume that the polynomi-als a and b are obtained by the shell reduction in (5). So thedecomposition (5) holds for the present notation. Moreover,let M K be the subspace of polynomial reduction w.r.t. K ,and N K its standard complement.We are going to determine necessary and sufficient condi-tions on hyperexponential integrability. Since gcd( b, k )=1, abk = p + qb + rk , (7)where p, q, r ∈ E [ y ], deg( q ) < deg( b ), and deg( r ) < deg( k ).Using an echelon basis of M K , we compute u in M K and v in N K s.t. k p + r = u + v. By the definition of M K , thereexists w in E [ y ] s.t. u = k D y ( w ) + k w . By (7), we get abk = qb + k D y ( w ) + k w + vk = D y ( w ) + Kw + qb + vk . It follows that abk T = D y ( wT ) + (cid:18) qb + vk (cid:19) T. (8)The process for obtaining (8) is referred to as the polynomialreduction for ( a/ ( bk )) T w.r.t. K , as it makes essential useof the subspaces M K and N K . By (8) and (5), H = D y (( S + w ) T ) + (cid:18) qb + vk (cid:19) T, (9)which motivates us to introduce the notion of residual forms. Definition 10.
With Convention 3, we further let f be arational function in E ( y ) . Another rational function r ∈ E ( y ) is said to be a residual form of f w.r.t. K if there exist g in E ( y ) and q, b, p in E [ y ] s.t. f = D y ( g ) + Kg + r and r = qb + vk , where b is squarefree, gcd( b, k ) = 1 , deg q < deg b , and v is in the standard complement N K of the subspace of poly-nomial reduction w.r.t. K . For brevity, we say that r is aresidual form w.r.t K if f is clear from context. Residual forms are closely related to the integrability ofhyperexponential functions.
Lemma 11.
With Convention 3, we further assume that r is a nonzero residual form w.r.t. K . Then the hyperexpo-nential function rT is not integrable.Proof. Let M K be the subspace for polynomial reduction,and N K its standard complement w.r.t. K . By the defini-tion of residual forms, there exist b, q ∈ E [ y ] with b beingsquarefree and v ∈ N K s.t.deg b > deg q, gcd( b, k ) = 1 , and r = qb + vk . (10)Thus, r can be rewritten as a/ ( bk ) for some a ∈ E [ y ]. Notethat a is not necessarily coprime with bk . It follows that rT = ab exp (cid:18)Z k − D y ( k ) k dy (cid:19) . Since ( k − D y ( k )) /k is differential-reduced and k , b arecoprime, ( a/b, ( k − D y ( k )) /k ) is indecomposable accord-ing to Definition 2 in [10]. By Theorem 4 in [10], a/b isin E [ y ]. So the denominator of r divides k , which, togetherwith (10), implies that q = 0. Consequently, ( v/k ) T is in-tegrable. By (2), v = k D y ( s ) + k s for some s ∈ E ( y ).Since v ∈ E [ y ], v ∈ M K by Lemma 6 (i). Thus, v = 0 be-cause v ∈ N K . We have that r = 0, a contradiction to theassumption that r = 0.The existence and uniqueness of residual forms are de-scribed below. Lemma 12.
With Convention 3, we have that the shell S has a residual form w.r.t. the kernel K . If a rational functionhas two residual forms w.r.t. K , then they are equal.Proof. By (9), S = D y ( S + w ) + ( S + w ) K + q/b + v/k . So q/b + v/k is a required form.Let r = q/b + v/k and r ′ = q ′ /b ′ + v ′ /k be two residualforms of a rational function w.r.t. K , where b, b ′ , q, q ′ , v, v ′ are in E [ y ], b and b ′ are squarefree, gcd( b, k )= gcd( b ′ , k )=1,eg q < deg b , deg q ′ < deg b ′ and v, v ′ ∈ N K . By the defini-tion of residual forms, D y ( f ) + fK + r = D y ( f ′ ) + f ′ K + r ′ for some f, f ′ ∈ E ( y ). It follows that D y (cid:0) f − f ′ (cid:1) + (cid:0) f − f ′ (cid:1) K + r − r ′ = 0 . Hence, ( r ′ − r ) T is integrable by (2). Since r − r ′ is also aresidual form w.r.t. K , r = r ′ by Lemma 11.Below is the main result of the present section. Theorem 13.
Let H be a hyperexponential function whoselogarithmic derivative has kernel K and shell S . Then thereis an algorithm for computing a rational function h in E ( y ) and a unique residual form r w.r.t. K s.t. H = D y (cid:18) h exp (cid:18)Z K dy (cid:19)(cid:19) + r exp (cid:18)Z K dy (cid:19) . (11) Moreover, H is integrable if and only if r = 0 .Proof. Let T = exp (cid:0)R K dy (cid:1) . Applying the shell reduc-tion to H w.r.t. K , we can find a rational function S ,and two polynomials a, b s.t. (4) holds. Then we applythe polynomial reduction to a/ ( bk ) T to get the residualform r = q/b + v/k s.t. (11) holds.Suppose that there exists another decomposition H = D y (cid:0) h ′ T (cid:1) + r ′ T (12)for some h ′ ∈ E ( y ) and r ′ is a residual form w.r.t. K . Thenboth r and r ′ are residual forms of S by (11), (12) and thefact H = ST . So r = r ′ by Lemma 12.If r = 0, then H is obviously integrable. Conversely,assume that H is integrable. Then rT is also integrableby (11). So r = 0 by Lemma 11.The reduction algorithm described in the proof of Theo-rem 13 has three interesting features. First, it enables us todecide hyperexponential integrability immediately. Second,it decomposes a hyperexponential function into a sum of anintegrable one and a non-integrable one in a canonical way.Third, it does not need to compute a polynomial solutionof any first-order linear differential equation. The methodwill be referred to as Hermite reduction for hyperexponen-tial functions in the sequel, because it extends all importantconclusions obtained by Hermite reduction for rational func-tions to hyperexponential ones.
Example 14.
Let H be the same hyperexponential functionas in Example 2. Then K = y/ ( y + 1) and S = 1 / ( y − .Set T = p y + 1 . By the shell reduction in Example 4, H = D y (cid:18) − y − T (cid:19) + ybk T, where b = y − and k = y + 1 . The polynomial reduc-tion yields ( y/ ( bk )) T = D y ( − T /
2) + (1 / (2 b ) + 1 / (2 k )) T. Combining the above equations, we decompose H as H = D y (cid:18) − ( y + 1)2( y − T (cid:19) + (cid:18) b + 12 k (cid:19) T. Example 15.
Consider H = y exp( y ) as given in Exam-ple 5. Since its logarithmic derivative has kernel K = 1 , thesubspace M K of polynomial reduction is equal to E [ y ] . Thus, y ∈ M K and H is integrable. More generally, M K = E [ y ] corresponds to the wellknown fact that p ( y ) exp( y ) is inte-grable for all p ∈ E [ y ] \ { } .
5. KERNEL REDUCTION
Let K = k /k be a nonzero differential-reduced rationalfunction in E ( y ) with gcd( k , k ) = 1. We may want toreduce a hyperexponential function in the form pk m exp (cid:18)Z K dy (cid:19) for some p ∈ E [ y ] and m ∈ N . One way would be to rewrite the above function as p exp (cid:18)Z k − mD y ( k ) k dy (cid:19) , and proceed by polynomial reduction w.r.t. the new ker-nel ( k − mD y ( k )) /k , which is also differential-reduced.However, it will prove to be more convenient in Section 6 toreduce the given function w.r.t. the initial kernel K . To thisend, we introduce another type of reduction, based on theideas in [9, 18]. Lemma 16.
With Convention 3, we let p ∈ E [ y ] and m ≥ .Then there exist p , p ∈ E [ y ] s.t. pk m = D y (cid:18) p k m − (cid:19) + p k m − K + p k . (13) Proof.
We proceed by induction on m . If m = 1, then tak-ing p = 0 and p = p yields the claimed form. Assumethat m >
1. We first show that there exist ˜ p , ˜ p ∈ E [ y ] s.t. pk m = D y (cid:18) ˜ p k m − (cid:19) + ˜ p k m − K + ˜ p k m − , which is equivalent to p = ˜ p ( k − ( m − D y ( k )) + ( D y (˜ p ) + ˜ p ) k . Since k /k is differential-reduced, there exist u, v ∈ E [ y ]s.t. p = u ( k − ( m − D y ( k )) + vk by the extended Eu-clidean algorithm. So we can take ˜ p = u and ˜ p = v − D y ( u ).By the induction hypothesis, there exist ¯ p , ¯ p ∈ E [ y ] s.t.˜ p k m − = D y (cid:18) ¯ p k m − (cid:19) + ¯ p k m − K + ¯ p k . Setting p = ¯ p k + ˜ p and p = ¯ p completes the proof.With Convention 3, we have pk m T = D y (cid:18) p k m − T (cid:19) + p k T by Lemma 16. This reduction will be referred to as the kernel reduction for ( p/k m ) T w.r.t. K .
6. TELESCOPING VIA REDUCTIONS
Hermite reduction has been used to construct telescopersfor bivariate rational functions in [5]. The goal of this sec-tion is to develop a reduction-based telescoping method forbivariate hyperexponential functions.
We briefly recall the reduction-based method for rational-function telescoping in [5].Let F be a field of characteristic zero and F ( x, y ) be thefield of rational functions in x and y over F . Let D x and D y denote the usual derivations ∂/∂x and ∂/∂y , respectively.et F ( x ) h D x i be the ring of linear differential operatorsover F ( x ). The ring F ( x ) h D x i is a left Euclidean domainand its left ideals are principal. For r ∈ F ( x, y ), the tele-scoping problem is to construct a nonzero linear differen-tial operator L ( x, D x ) ∈ F ( x ) h D x i s.t. L ( x, D x )( r ) = D y ( s ) , where s ∈ F ( x, y ). The operator L is called a telescoper for r , and s is the corresponding certificate . The set T ofall telescopers for a given rational function is a left idealof F ( x ) h D x i . Any generator of T is called a minimal tele-scoper for the given rational function.For any i ∈ N , rational Hermite reduction (w.r.t. y ) de-composes D ix ( r ) into D ix ( r )= D y ( s i )+ a i /b, where s i ∈ F ( x, y )and a i , b ∈ F ( x )[ y ] with deg y ( a i ) < deg y ( b ), and b is square-free over F ( x ). Since deg y ( a i ) is bounded by deg y ( b ), the se-quence { a i } i ∈ N is linearly dependent over F ( x ). Assume thatthere exist e , . . . , e ρ ∈ F ( x ), not all zero, s.t. P ρi =0 e i a i = 0.Then L := P ρi =0 e i D ix is a telescoper for r and P ρi =0 e i g i isthe corresponding certificate. In fact, L is a minimal tele-scoper for r if ρ is the minimal integer s.t. e , . . . , e ρ arelinearly dependent over F ( x ). This reasoning yields the up-per bound deg y ( b ) on the order of minimal telescopers. We now apply the Hermite reduction for univariate hyper-exponential functions in Section 4 to compute telescopers forbivariate hyperexponential functions.A nonzero element H in some differential field extensionof F ( x, y ) is said to be hyperexponential over F ( x, y ) if its log-arithmic derivatives D x ( H ) /H and D y ( H ) /H are in F ( x, y ).Put f = D x ( H ) /H and g = D y ( H ) /H . Then D y ( f )= D x ( g )because D x and D y commute. Therefore, it is legitimate todenote H by exp( R f dx + g dy ). For two hyperexponentialfunctions H i = exp( R f i dx + g i dy ), i = 1 ,
2, we have H H = exp (cid:18)Z ( f + f ) dx + ( g + g ) dy (cid:19) . In particular, the product of a rational function r ∈ F ( x, y )and a hyperexponential function H = exp( R f dx + g dy ) is rH = exp (cid:18)Z ( f + D x ( r ) /r ) dx + ( g + D y ( r ) /r ) dy (cid:19) . The following fact is immediate from [10, Lemma 8].
Fact 17.
Let f and g be rational functions in F ( x, y ) sat-isfying D y ( f ) = D x ( g ) . Then the denominator of f dividesthat of g in F ( x )[ y ] . For a hyperexponential function H over F ( x, y ), the tele-scoping problem is to construct a linear ordinary differentialoperator L ( x, D x ) in F ( x ) h D x i s.t. L ( x, D x )( H ) = D y ( G )for some hyperexponential function G over F ( x, y ). As in therational case, our idea is to apply the Hermite reduction forunivariate hyperexponential functions w.r.t. y to the deriva-tives D ix ( H ) iteratively, and then find a linear dependencyamong the residual forms over F ( x ). Lemma 18.
Let H = exp( R f dx + g dy ) be a hyperexpo-nential function over F ( x, y ) . Let K be the kernel and S theshell of g w.r.t. y . Then, for every i ∈ N , the i -th deriva-tive D ix ( H ) can be decomposed into D ix ( H ) = D y ( u i T ) + r i T, (14) where u i ∈ F ( x, y ) , T = exp( R ( f − D x ( S ) /S ) dx + K dy ) and r i ∈ F ( x, y ) is a residual form w.r.t. K . Moreover,let k be the denominator of K , b the squarefree part of thedenominator of S , and N K the standard complement of thesubspace for polynomial reduction w.r.t. K . Then r i = q i b + v i k (15) for some q i ∈ F ( x )[ y ] with deg y q i < deg y b and v i ∈ N K .Proof. We proceed by induction on i . If i = 0, then theassertion holds by Theorem 13.Assume that D ix ( H ) can be decomposed into (14) andassume that (15) holds. Moreover, let ˜ f = f − D x ( S ) /S .Consider the ( i + 1)-th derivative D i +1 x ( h ). There exists apolynomial a in F ( x )[ y ] s.t. ˜ f = a/k by D y (cid:16) ˜ f (cid:17) = D x ( K )and Fact 17. A direct calculation leads to D i +1 x ( H ) = D y ( D x ( u i T )) + (cid:18) aq i bk + D x ( q i ) b + D x ( v i ) k (cid:19) T + (cid:18) − q i D x ( b ) b + ( a − D x ( k )) v i k (cid:19) T. Applying the shell reduction to (cid:0) − q i D x ( b ) /b (cid:1) T and thekernel reduction to (cid:0) ( a − D x ( k )) v i /k (cid:1) T w.r.t. y , we get − q i D x ( b ) b = D y (cid:16) w b (cid:17) + w b K + w bk , ( a − D x ( k )) v i k = D y (cid:18) p k (cid:19) + p k K + p k , where w , w , p and p are in F ( x )[ y ]. We then apply poly-nomial reduction to ˜ ST w.r.t. K , where˜ S = w bk + p k + aq i bk + D x ( q i ) b + D x ( v i ) k , which leads to˜ S = D y ( w ) + wK + (cid:18) q i +1 b + v i +1 k (cid:19) , where w ∈ F ( x, y ) and q i +1 /b + v i +1 /k is the residual formof ˜ S w.r.t. K . It follows from a direct calculation that D i +1 x ( H ) = D y ( u i +1 T ) + (cid:18) q i +1 b + v i +1 k (cid:19) T, where u i +1 = D x ( u i ) + u i ˜ f + w /b + p /k + w .The main results in the present section are given below. Theorem 19.
With the notation introduced in Lemma 18,we let L = P ρi =0 e i D iy with e , . . . , e ρ ∈ F ( x ) , not all zero.(i) L is a telescoper for H if and only if P ρi =1 e i r i = 0 . (ii) The order of a minimal telescoper for H is no morethan deg y ( b ) + max(deg y ( k ) , deg y ( k ) − .Proof. We set E = F ( x ) and view that hyperexponentialfunctions involved in the proof are over E ( y ). Moreover,let u = P ρi =0 e i u i and r = P ρi =0 e i r i . By (14), we have L ( H ) = D y ( uT ) + rT. (16)If r = 0, then L is a telescoper by (16). Conversely, assumethat L is a telescoper of h . Then rT is integrable w.r.t. y y (16). Since r is a residual form, it is equal to zero byLemma 11. The first assertion is proved.Set λ = max(deg y ( k ) , deg y ( k ) − r i = q i /b + v i /k be as defined in (14) and (15). ByCorollary 9, the v i ’s have a common set P of supportingmonomials with |P| ≤ λ . Moreover, deg y ( q i ) < deg y ( b )and gcd( b, k ) = 1. Therefore, the residual forms r , . . . , r ρ are linearly dependent over F ( x ) if ρ ≥ deg y ( b ) + λ . Thesecond assertion holds Remark 20.
By Theorem 19, the first linear dependencyamong the residual forms r , r , r , . . . gives rise to a mini-mal telescoper of H . Below is an outline of the reduction based telescoping algo-rithm for hyperexponential functions, in which the notationis that introduced in Lemma 18 is used.
Algorithm.
HermiteTelescoping : Given a bivariate hyper-exponential function H = exp( R f dx + g dy ) over F ( x, y ),compute a minimal telescoper L and its certificate w.r.t. y .1. Find the kernel K and shell S of D y ( H ) /H w.r.t. y .Set b to be the squarefree part of the denominator of S .2. Decompose H into H = D y ( u T ) + r T using the Her-mite reduction for hyperexponential functions given inTheorem 13. If r = 0, return (1 , u T ).3. Set ρ := deg y ( b ) + max(deg y ( k ) , deg y ( k ) − i from 0 to ρ do4.1. Compute ( u i , r i ) incrementally s.t. D ix ( H ) = D y ( u i T ) + r i T by the shell, kernel and polynomial reductions de-scribed in Lemma 18.4.2. Find η j ∈ F ( x ) s.t. P ij =0 η j r j = 0 using the al-gorithm in [17]. If there is a nontrivial solution,return (cid:16)P ij =0 η j D jx , P ij =0 η j u j T (cid:17) . Example 21.
Let H = √ x − y exp( x y ) . Then D x ( H ) /H and D y ( H ) /H are, respectively, f = 1 + 4 x y − xy x − y ) and g = − x − x yx − y . Since g is differential-reduced w.r.t. y , g is the kernel and is the shell of D y ( H ) /H w.r.t. y . By Hermite reduction, H = D y (cid:18) x H (cid:19) + 1 x k H. (17) Applying D x to the above equation yields D x ( H ) = D y (cid:18) − x + 8 y + 4 x y − x y x ( x − y ) H (cid:19) + rH, where r = ( − x + 8 y + 4 x y − x y ) / (2 x k ) . The shell,kernel and polynomial reduction given in Lemma 18 yields D x ( H ) = D y (cid:18) x y − x · H (cid:19) + 3 x − x k H (18) Combining (17) and (18) , we get L = (6 − x ) + 2 xD x is a minimal telescoper for H and G = (4 y − x ) H is thecorresponding certificate. Remark 22.
The algorithm
HermiteTelescoping is directlybased on the proof of Lemma 18. Yet, there is another ideafor computing a minimal telescoper of H . Namely, we firstcompute a nonzero operator L ∈ F ( x ) h D x i of minimal orders.t. L ( H ) = D y ( G ) + ( p/k ) T for some hyperexponentialfunction G and polynomial p . Note that such operators al-ways exist, because deg y q i in (15) is less than deg y b . Thenwe apply the algorithm HermiteTelescoping to get a minimaltelescoper L for ( p/k ) T . In doing so, any rational functionwith denominator b will not appear when we compute L . Itturns out that L L is a minimal telescoper of H . An im-plementation on this idea is underway. Assume that H = u exp (cid:18) r r (cid:19) m Y i =1 p i ( x, y ) c i , (19)where u, r , r , p , . . . , p m are nonzero polynomials in F [ x, y ]and c , . . . , c m are distinct indeterminates . Theorem cAZin [4] asserts that the order of minimal telescopers for H isbounded by α := deg y ( r ) + max (cid:0) deg y ( r ) , deg y ( r ) (cid:1) + m X i =1 deg y ( p i ) − . Note that H can be viewed as a hyperexponential functionover F ( c , . . . , c m )( x, y ). We now show that α given above isno less than the order bound on minimal telescopers for H obtained from Theorem 19 (ii). The kernel and shell of thelogarithmic derivative D y ( H ) /H are K := D y (cid:18) r r (cid:19) + m X i =1 c i D y ( p i ) p i and S := u, respectively, because K has no integral residue at any simplepole, S is a polynomial in F [ x, y ], and D y ( H ) /H is equalto K + D y ( S ) /S . Let K = k /k with gcd( k , k ) = 1. Adirect calculation leads todeg y ( k ) ≤ deg y ( r ) + deg y ( r ) + m X i =1 deg y ( p i ) − , and deg y ( k ) ≤ y ( r ) + m X i =1 deg y ( p i ) . By Theorem 19, the order of minimal telescopers for H is nomore than max (cid:0) deg y ( k ) , deg y ( k ) − (cid:1) , which is no morethan α by the above two inequalities.Indeed, the order bound in Theorem 19 (ii) may be smallerthan that in Theorem cAZ. Example 23.
Let H = q c exp( a/q ) , where a, q are irreduciblepolynomials in F [ x, y ] with deg y ( a ) < deg y ( q ) , and c is atranscendental constant over F . By Theorem cAZ, a min-imal telescoper for H has order no more than y q − .On the other hand, the kernel and shell of D y ( H ) /H areequal to ( D y ( a ) q − aD y ( q )+ cqD y ( q )) /q and , respectively.A minimal telescoper has order no more than y q − byTheorem 19 (ii). n general, Christopher’s Theorem states that a hyperex-ponential function over F ( x, y ) can always be written as: uv exp (cid:18) r r (cid:19) m Y i =1 p i ( x, y ) c i , (20)where u, v, r , r ∈ F [ x, y ], c i is algebraic over F , and p i is in F ( c i )[ x, y ], i = 1 , . . . , m . A more explicit descriptionon (20) can be found in [7]. So H given in (19) is a spe-cial instance for hyperexponential functions. In addition, itis easier to compute the kernel and shell than to computethe decompositions (19) and (20) when a hyperexponentialfunction is given by its logarithmic derivatives.
7. IMPLEMENTATION AND TIMINGS
We have produced a preliminary implementation of thealgorithm
HermiteTelescoping in the computer algebra sys-tem
Maple 16 . Our Maple code is available from
We now compare the performance of our algorithm to theMaple function
DEtools[Zeilberger] of the telescoping algo-rithm in [3]. The examples for comparison are of the form pq m · r ab · exp (cid:16) uv (cid:17) , where m ∈ N , p, q, a, b, u, v ∈ Z [ x, y ] are irreducible andtheir coefficients are randomly chosen. For simplicity, wechoose λ = deg y ( p ) = deg y ( q ), µ = deg y ( a ) = deg y ( b ),and ν = deg y ( u ) = deg y ( v ). The runtime comparison (inseconds) for different examples is shown in Table 1, in which • ZT : the Maple function DEtools[Zeilberger] . • HT : our implementation of HermiteTelescoping . • order: the order of the computed minimal telescoper. • OOM: Maple runs out of memory.( λ, µ, ν, m ) ZT HT order(2, 0, 2, 1) 2.23 2.43 5(2, 0, 2, 2) 2.21 2.01 5(3, 0, 2, 1) 8.72 6.64 6(3, 0, 2, 2) 9.38 6.56 6(6, 0, 1, 1) 45.35 24.49 7(6, 0, 1, 2) 43.02 22.91 7(2, 2, 2, 1) 1405.9 221.5 9(2, 2, 2, 2) 1398.1 200.34 9(3, 0, 3, 1) 147.92 47.23 8(3, 0, 3, 2) 151.20 44.56 8(3, 3, 0, 1) 207.82 61.10 8(3, 3, 0, 2) 211.70 58.63 8(3, 2, 1, 1) 304.61 62.67 8(3, 2, 1, 2) 331.0 63.61 8(3, 1, 3, 1) OOM 534.87 10(3, 1, 3, 2) OOM 522.15 10
Table 1:
Timings (in sec.) were taken on a Mac OSX computer with 4Gb RAM and 3.06 GHz Core 2 Duoprocessor.
Remark 24.
The orders of the computed minimal tele-scopers in our experiments are equal to the predicted orderbounds in Theorem 19 .
8. REFERENCES [1] S.A. Abramov. The rational component of the solution of afirst order linear recurrence relation with rational righthand-side. ˇZ. Vyˇcisl. Mat. i Mat. Fiz. , 15(4):1035–1039,1090, 1975.[2] S.A. Abramov and M. Petkovˇsek. Minimal decompositionof indefinite hypergeometric sums. In
ISSAC’01:Proceedings of the 2001 International Symposium onSymbolic and Algebraic Computation , pages 7–14, NewYork, 2001. ACM.[3] G. Almkvist and D. Zeilberger. The method ofdifferentiating under the integral sign.
J. SymbolicComput. , 10:571–591, 1990.[4] M. Apagodu and D. Zeilberger. Multi-variable Zeilbergerand Almkvist-Zeilberger algorithms and the sharpening ofWilf- Zeilberger theory.
Adv. in Appl. Math. ,37(2):139–152, 2006.[5] A. Bostan, S. Chen, F. Chyzak, and Z. Li. Complexity ofcreative telescoping for bivariate rational functions. In
ISSAC’10: Proceedings of the 2010 InternationalSymposium on Symbolic and Algebraic Computation , pages203–210, New York, NY, USA, 2010. ACM.[6] M. Bronstein.
Symbolic Integration I: TranscendentalFunctions , volume 1 of
Algorithms and Computation inMathematics . Springer-Verlag, Berlin, second edition, 2005.[7] S. Chen.
Some applications of differential-differencealgebra to creative telescoping . PhD thesis, ´EcolePolytechnique (Palaiseau, France), February 2011.[8] C. Christopher. Liouvillian first integrals of second orderpolynomial differential equations.
Electron. J. DifferentialEquations , 49:1–7, 1999.[9] J.H. Davenport. The Risch differential equation problem.
SIAM J. Comput. , 15(4):903–918, 1986.[10] K.O. Geddes, H.Q. Le, and Z. Li. Differential rationalnormal forms and a reduction algorithm forhyperexponential functions. In
ISSAC’04: Proceedings ofthe 2004 International Symposium on Symbolic andAlgebraic Computation , pages 183–190, New York, USA,2004. ACM.[11] R.W. Gosper, Jr. Decision procedure for indefinitehypergeometric summation.
Proc. Nat. Acad. Sci. U.S.A. ,75(1):40–42, 1978.[12] C. Hermite. Sur l’int´egration des fractions rationnelles.
Ann. Sci. ´Ecole Norm. Sup. (2) , 1:215–218, 1872.[13] E. Horowitz. Algorithms for partial fraction decompositionand rational function integration. In
SYMSAC’71 , pages441–457, New York, USA, 1971. ACM.[14] M. Kauers, C. Koutschan, and D. Zeilberger. Proof of IraGessel’s lattice path conjecture.
Proc. Natl. Acad. Sci.USA , 106(28):11502–11505, 2009.[15] C. Koutschan, M. Kauers, and D. Zeilberger. Proof ofGeorge Andrews’s and David Robbins’s q -TSPP conjecture. Proc. Natl. Acad. Sci. USA , 108(6):2196–2199, 2011.[16] M.V. Ostrogradski˘ı. De l’int´egration des fractionsrationnelles.
Bull. de la classe physico-math´ematique del’Acad. Imp´eriale des Sciences de Saint-P´etersbourg ,4:145–167, 286–300, 1845.[17] A. Storjohann and G. Villard. Computing the rank and asmall nullspace basis of a polynomial matrix. In
ISSAC’05:Proceedings of the 2005 International Symposium onSymbolic and Algebraic Computation , pages 309–316.ACM, New York, 2005.[18] G. Xin and T.Y.J. Zhang. Enumeration of bilaterallysymmetric 3-noncrossing partitions.