Efficient q-Integer Linear Decomposition of Multivariate Polynomials
aa r X i v : . [ c s . S C ] F e b E ffi cient q -Integer Linear Decompositionof Multivariate Polynomials Mark Giesbrecht
Symbolic Computation Group, Cheriton School of Computer Science, University of Waterloo,Waterloo, ON, N2L 3G1, Canada
Hui Huang
Symbolic Computation Group, Cheriton School of Computer Science, University of Waterloo,Waterloo, ON, N2L 3G1, Canada
George Labahn
Symbolic Computation Group, Cheriton School of Computer Science, University of Waterloo,Waterloo, ON, N2L 3G1, Canada
Eugene Zima
Physics and Computer Science, Wilfrid Laurier University,Waterloo, ON, N2L 3C5, Canada
Abstract
We present two new algorithms for the computation of the q -integer linear decomposition of amultivariate polynomial. Such a decomposition is essential in the q -analogous world of symbolicsummation, for example, describing the q -counterpart of Ore-Sato theory or determining the ap-plicability of the q -analogue of Zeilberger’s algorithm to a q -hypergeometric term. Both of ouralgorithms require only basic integer and polynomial arithmetic and work for any unique factor-ization domain containing the ring of integers. Complete complexity analyses are conducted forboth our algorithms and two previous algorithms in the case of multivariate integer polynomials,showing that our algorithms have better theoretical performances. A Maple implementation isalso included which suggests that our algorithms are also much faster in practice than previousalgorithms. Keywords: q -Analogue, Integer-linear polynomials, Polynomial decomposition,Newton polytope (Newton polygon), Creative telescoping, Ore-Sato theory
1. Introduction
Many objects in the ordinary shift world of symbolic summation find a natural counterpartcommonly called q-analogues . In a typical situation, these are just slight adaptations of the
Email addresses: [email protected] (Mark Giesbrecht), [email protected] (Hui Huang), [email protected] (George Labahn), [email protected] (Eugene Zima) riginal objects but with involved variables promoted to exponents of an additional parameter q . Techniques for handling the originals often carry over to their q -analogues with some subtlemodifications.In this paper, we deal with the q -analogue of integer-linear decompositions of polynomialsand aim to provide an intensive treatment for its computation in analogy to (Giesbrecht et al.,2019a). Surprisingly, although this q -analogue is obtained by modeling its ordinary shift coun-terpart, the primary technique used in (Giesbrecht et al., 2019a) can not be easily adapted tocompute it due to di ff erent structures. A new alternative technique will be presented in this q -shift case.In order to describe more details, we let D be a ring of characteristic zero and let R = D [ q , q − ] be its transcendental ring extension by the indeterminate q . For n discrete indeter-minates k , . . . , k n distinct from q , we know that q k , . . . , q k n are transcendental over R . We canthen consider polynomials in q k , . . . , q k n over R , all of which form a well-defined ring denotedby R [ q k , . . . , q k n ]. We say an irreducible polynomial p ∈ R [ q k , . . . , q k n ] is q-integer linear over R if there exists a univariate polynomial P ( y ) ∈ R [ y ] and two integer-linear polynomials P ni = α i k i , P ni = λ i k i ∈ Z [ k , . . . , k n ] such that p ( q k , . . . , q k n ) = q P ni = α i k i P ( q P ni = λ i k i ) . In order to avoid superscripts, we will write the indeterminates q k , . . . , q k n as the variables x , . . . , x n in the sequel of the paper. Then the above definition can be rephrased as follows.An irreducible polynomial p ∈ R [ x , . . . , x n ] is called q-integer linear over R if there exists aunivariate polynomial P ( y ) ∈ R [ y ] and integers α , . . . , α n , λ , . . . , λ n such that p ( x , . . . , x n ) = x α · · · x α n n P ( x λ · · · x λ n n ) . (1.1)Note that the indeterminate q is hidden in the variables x , . . . , x n . Since a common factor ofthe λ i can be pulled out and absorbed into P , and a monomial can be merged into x α · · · x α n n ifnecessary, we assume that the integers λ , . . . , λ n have no common divisor, that the last nonzerointeger is nonnegative, that λ i = x i ( p ) = P (0) ,
0. Such a vector( λ , . . . , λ n ), as well as such a polynomial P ( y ), is unique. We call the vector ( λ , . . . , λ n ) the q-integer linear type of p and the polynomial P ( y ) its corresponding univariate polynomial .Note that the resulting α , . . . , α n all belong to N since p ∈ R [ x , . . . , x n ] and P ( y ) ∈ R [ y ]. Apolynomial in R [ x , . . . , x n ] is called q-integer linear (over R ) if all its irreducible factors are q -integer linear, possibly with di ff erent q -integer linear types. For a polynomial p ∈ R [ x , . . . , x n ],we can define its q-integer linear decomposition by factoring into irreducible q -integer linear ornon- q -integer linear polynomials and collecting irreducible factors having common types.As with its ordinary shift counterpart, q -integer linear polynomials find broad applicationsin the q -analysis of symbolic summation. In particular, it is an important ingredient of the q -analogue of the Ore-Sato theorem for describing the structure of multivariate q -hypergeometricterms (Du and Li, 2019), which in turn, as indicated by (Chen and Koutschan, 2019), serves asa promising indispensable tool for settling a q -analogue of Wilf-Zeilberger’s conjecture givenin (Wilf and Zeilberger, 1992). Moreover, the q -integer linearity of polynomials plays a crucialrole in detecting the applicability of the q -analogue of Zeilberger’s algorithm (also known as themethod of creative telescoping) for bivariate q -hypergeometric terms (Chen et al., 2005).The full q -integer linear decomposition of polynomials is also very useful. On the one hand,it provides a natural way to determine the q -integer linearity of a given polynomial. On the2ther hand, it enables one to compute the q -analogue of Ore-Sato decomposition of a given q -hypergeometic term, and can also be employed to develop a fast creative telescoping algorithmfor rational functions in the q -shift setting in analogy to (Giesbrecht et al., 2019b). Evidently, thee ffi ciency of the computation of q -integer linear decompositions directly a ff ects the utility of allthese algorithms.In contrast to the ordinary shift case (cf. (Abramov and Le, 2002; Giesbrecht et al., 2019a;Li and Zhang, 2013)), algorithms for computing the q -integer linear decomposition of a mul-tivariate polynomial are not very well developed. As far as we are aware, there is only onealgorithm available to compute such a decomposition of a bivariate polynomial. This algorithmwas first developed by Le (2001, §
5) with an extended description provided in (Le et al., 2001).Except for using the same pattern as its ordinary shift counterpart (Abramov and Le, 2002), thisalgorithm takes use of a completely di ff erent strategy, especially for finding q -integer lineartypes. This is mainly because all q -integer linear types appear as the exponent vectors of p ,rather than as the coe ffi cients in the ordinary shift case. The main idea used by Le (2001, §
5) isto first find candidates for q -integer linear types by computing a resultant and then, for each can-didate, extract the corresponding univariate polynomial via bivariate GCD computations. Giventhe algebraic machinery on which the algorithm is based, it is not clear how one can directlygeneralize this to handle polynomials in more than two variables.Our first contribution is a bivariate-based scheme, similar to the one given in (Giesbrecht et al.,2019a, § q -analogue of the algorithm developed by Li and Zhang(2013). This algorithm makes use of the observation that the di ff erence of exponent vectors ofany two monomials appearing in an irreducible q -integer linear polynomial, say the polynomial p of the form (1.1), must be a scalar multiple of the q -integer linear type ( λ , . . . , λ n ). We alsogive a complexity analysis for both algorithms, at least in the case of bivariate polynomials over Z [ q , q − ], supporting the superiority of the factorization-based algorithm over the algorithm ofLe. Same opinion is suggested by the empirical tests.The main contribution of this paper is a pair of new fast algorithms for computing the q -integer linear decomposition of a multivariate polynomial. Both algorithms will work for anyunique factorization domain containing all integers and for any polynomial with an arbitrarynumber of variables. The first approach combines the main ideas of the two previous algorithms- in the sense that it follows the pattern of the algorithm of Le and also makes use of exponentvectors as the factorization-based algorithm - but avoids the computation of resultants as well asthe need for full irreducible factorization. More precisely, this approach reduces the problem offinding candidates for q -integer linear types to an easy ad hoc geometry task of constructing theNewton polytope of the given polynomial, implying computations only using basic arithmeticoperations ( + , − , ÷ , × ) of integers, and then computes each corresponding univariate polynomialby a simple content computation as in the ordinary shift case. As such we show that the q -analogue is actually simpler than its ordinary shift counterpart in the sense that, instead of findingrational roots of polynomials, one merely needs to perform basic integer manipulations.Our second approach relies on the bivariate-based scheme mentioned earlier. This schemetakes any bivariate algorithm, that is, an algorithm for computing the q -integer linear decom-position of a bivariate polynomial, as a base case and iteratively tackles only two variables at atime until all variables are treated. This scheme, together with the bivariate restriction of the firstapproach, immediately establishes our second approach. Clearly, these two approaches coincidein the bivariate case. 3oth approaches appear to be e ffi cient in practice, though the second method shows an ad-vantage for polynomials of a large number of variables. In order to do a theoretical comparisonwe have analyzed the worst-case running time complexity of both approaches in the case of mul-tivariate polynomials over Z [ q , q − ]. The analysis shows that the second approach is superior tothe first one when the given polynomial has more than two variables. When restricted to the caseof bivariate polynomials over Z [ q , q − ], the two approaches merge into one, which in turn is con-siderably faster than the algorithm of Le and the algorithm based on factorization. In addition,we also give experimental results which verify our complexity comparisons.The remainder of the paper proceeds as follows. Background and basic notions required inthe paper are provided in the next section. In Section 3 we present a fast algorithm for computingthe q -integer linear decomposition in the special case of a bivariate polynomial including itscomplexity costs. Two approaches for the general multivariate case, along with their respectivecomplexity analyses, are given successively in Sections 4 and 5. The following section providesa complexity comparison of our bivariate algorithm, the algorithm of Le and the factorization-based algorithm. The paper ends with an experimental comparison among all algorithms, alongwith a conclusion section.
2. Preliminaries
Throughout the paper, we let D be a unique factorization domain (UFD) of characteristiczero with R = D [ q , q − ] denoting the transcendental ring extension by an indeterminate q . Notethat a domain of characteristic zero always contains the ring of integers Z as a subdomain. Let R [ x , . . . , x n ] be the ring of polynomials in x , . . . , x n over R , where x , . . . , x n are variablesdistinct from q . We reserve the variables x and y as synonyms for x and x , respectively, so asto avoid subscripts in the case when n ≤ p be a polynomial in R [ x , . . . , x n ]. Throughout this paper we will order monomials in R [ x , . . . , x n ] using a pure lexicographic order in x ≺ · · · ≺ x n . For this order we let lc( p )and deg( p ) denote the leading coe ffi cient and the total degree, respectively, of p with respectto x , . . . , x n . We follow the convention that deg(0) = −∞ . We say that p is monic (over R ) iflc( p ) =
1. The content of p (over R ), denoted by cont( p ), is the greatest common divisor (GCD)over R of the coe ffi cients of p with respect to x , . . . , x n with p being primitive if cont( p ) = primitive part prim( p ) of p (over R ) is defined as p / cont( p ). For brevity, we will omit thedomain if it is clear from the context. In certain instances, we also need to consider the abovenotions with respect to a subset of the n variables. In these cases, we will either specify therelevant domain or indicate the related variables as subscripts of the corresponding notion. Forexample, lc x , x ( p ), deg x , x ( p ), cont x , x ( p ) and prim x , x ( p ) denote each function but applied toa polynomial p viewing it as a polynomial in x , x over the domain R [ x , . . . , x n ].In order to obtain a canonical representation, we introduce the notion of q -primitive poly-nomials in the univariate case. A polynomial p ∈ R [ x ] is called q-primitive if it is primitiveover R and its constant term p (0) is nonzero. Note that this concept is a ring counterpart of q -monic polynomials introduced by Paule and Riese (1997). Clearly, any factor of a q -primitivepolynomial in R [ x ] is again q -primitive.We work with n -dimensional vectors from the lattice Z n . In order to simplify notations, weemploy bold letters, say i , for an n -dimensional vector ( i , . . . , i n ) ∈ Z n , and the multi-indexconvention x i for the monomial x i · · · x i n n . The zero n -dimensional vector is denoted by boldface . As usual, the operations, especially addition / subtraction and scalar multiplication, on vectorsin the lattice Z n are performed componentwise. Notice that there is a one-to-one correspondence4etween monomials in R [ x , . . . , x n ] and vectors in N n . From the monomial ordering x ≺ · · · ≺ x n , we then establish a natural total ordering “ < ” on N n in such a way that for any two vectors i , j ∈ N n , if x i ≺ x j then we also say that i < j .The support of polynomials plays a crucial role in our algorithms. Recall that the support of apolynomial p ∈ R [ x , . . . , x n ], denoted by supp( p ), is defined as the set of indices i ∈ N n with theproperty that the coe ffi cient of x i in p is nonzero. Roughly speaking, supp( p ) records exponentvectors of all monomials present in p . Clearly, the support of a polynomial in R [ x , . . . , x n ] is afinite set in N n , and it is the empty set if and only if the given polynomial is zero. As describedpreviously, we sometimes only need to consider supports of polynomials with respect to partof the variables and these instances will be identified by explicitly pointing out the involvedvariables.We are interested in finding the following decomposition of a polynomial, something brieflyalluded to in the introduction. Definition 2.1.
Let p ∈ R [ x , . . . , x n ] be a polynomial admitting the decompositionp = c x α P m Y i = P i ( x λ i ) , (2.1) where c ∈ R , m ∈ N , α ∈ N n , λ i = ( λ i , . . . , λ in ) ∈ Z n , P ∈ R [ x , . . . , x n ] and P i ( y ) ∈ R [ y ] . Then (2.1) is called the q -integer linear decomposition of p (over R ) if(1) P is primitive and none of its nonconstant irreducible factors is q-integer linear;(2) each P i ( y ) is q-primitive and of positive degree in y;(3) each λ i is a nonzero q-integer linear type, in other words, gcd( λ i , . . . , λ in ) = and itsrightmost nonzero entry is positive;(4) any two vectors from the λ i are distinct.We say that the λ i are q -integer linear types of p and each P i ( y ) is the corresponding univariatepolynomial of the type λ i . Evidently, p is q -integer linear if and only if P is a unit of R in the decomposition (2.1).Hence all univariate polynomials are q -integer linear. By full factorization, we see that everypolynomial admits a q -integer linear decomposition. Moreover, this decomposition is uniqueup to the order of factors and multiplication by units of R , according to the uniqueness of the q -integer linear types and the full factorization.We will be considering Laurent polynomials in the ring R [ x , x − , . . . , x n − , x − n − , x n ], namelypolynomials of the form x − α · · · x − α n − n − · p for some α , . . . , α n − ∈ N and p ∈ R [ x , . . . , x n ].More generally, let K throughout denote the quotient field of R . Then we are interested inpolynomials in x n over the field K ( x , . . . , x n − ), all of which form the ring K ( x , . . . , x n − )[ x n ],containing R [ x , x − , . . . , x n − , x − n − , x n ] as a subring. It is convenient to extend the definitionof content and primitive part to polynomials in this setting. Let p ∈ K ( x , . . . , x n − )[ x n ] be ofthe form P di = ( a i / b ) x in for d ∈ N and a i , b ∈ R [ x , . . . , x n − ]. Then we let the content cont x n ( p )of p with respect to x n be gcd( a , . . . , a d ) / b and the corresponding primitive part prim x n ( p ) = p / cont x n ( p ). Evidently, prim x n ( p ) ∈ R [ x , . . . , x n ]. In particular, if p ∈ K [ x , . . . , x n ] then wecan further consider its content with respect to all variables x , . . . , x n . To be specific, by writing5 = P i ∈ N n ( a i ,..., i n / b ) x i for a i ,..., i n , b ∈ R , we let its content cont( p ) (over K ) be a / b with a being the GCD of the a i ,..., i n over R and the primitive part prim( p ) = p / cont( p ). Note that thedefinition of leading coe ffi cient and degree extends to polynomials in K [ x , . . . , x n ] in a naturalmanner. Again, as before, all these notions can be used in terms of partial variables.Two basic lemmas are given below for later use. Lemma 2.2.
Let P ( y ) ∈ R [ y ] \ R with P (0) , and let λ ∈ Z n with gcd( λ , . . . , λ n ) = , λ , . . . , λ n − not all zero and λ n > . Let f ∈ K ( x , . . . , x n − )[ x n ] be a factor of P ( x λ ) which ismonic with respect to x n . Then(i) there exists c ∈ K , α , . . . , α n − ∈ Z and a factor g ( y ) ∈ R [ y ] of P ( y ) such that f = cx α · · · x α n − n − g ( x λ ) . Moreover, if f is of positive degree in x n then deg( g ( y )) > .(ii) P ( y ) is irreducible over R if and only if P ( x λ ) is irreducible over K ( x , . . . , x n − ) if andonly if prim x n ( P ( x λ )) is irreducible over R .Proof. (i) Since P (0) ,
0, all its roots in the algebraic closure K of K are nonzero. In order toprove the assertion, it is su ffi cient to show that x λ − r for any root r ∈ K of P ( y ) is irreducible over K ( x , . . . , x n − ), because then since f ∈ K ( x , . . . , x n − )[ x n ] is a factor of P ( x λ ) and it is monicwith respect to x n , it factors completely into irreducibles in K ( x , . . . , x n − )[ x n ] as follows f = s Y i = ( x − λ · · · x − λ n − n − )( x λ − r i ) = ( x − λ · · · x − λ n − n − ) s s Y i = ( x λ − r i ) , where s ∈ N with s ≤ deg( P ( y )) and the r i ∈ K are roots of P ( y ), and thus the assertion directlyfollows by pulling out the content over K .Let r ∈ K be a root of P ( y ) and suppose that x λ − r is reducible over K ( x , . . . , x n − ). Considerthe algebraic closure K ( x , . . . , x n − ) of K ( x , . . . , x n − ) and let ω ∈ K be a λ n -th root of unitysuch that ω λ n =
1. Since r is nonzero, the complete factorization of x λ − r over K ( x , . . . , x n − )is given by x λ − r = x λ · · · x λ n − n − λ n − Y i = (cid:16) x n − ω i r /λ n x − λ /λ n · · · x − λ n − /λ n n − (cid:17) . It then follows from the reducibility of x λ − r over K ( x , . . . , x n − ) that there exist i , . . . , i k ∈{ , . . . , λ n − } with 0 < k < λ n such that k Y j = (cid:16) x n − ω i j r /λ n x λ /λ n · · · x λ n − /λ n n − (cid:17) ∈ K ( x , . . . , x n − )[ x n ] . This implies that ( λ i /λ n ) k ∈ Z for all i = , . . . , n −
1. Thus λ n divides k · gcd( λ , . . . , λ n − )in Z . Since λ , . . . , λ n − are not all zero and gcd( λ , . . . , λ n ) =
1, we have λ n divides k in Z , acontradiction with 0 < k < λ n .(ii) For the first equivalence, the su ffi ciency is evident. In order to show the necessity, supposethat P ( x λ ) is reducible over K ( x , . . . , x n − ). Then there exists a factor f ∈ K ( x , . . . , x n − )[ x n ]of P ( x λ ) which is of positive degree in x n . By the assertion (i), we then obtain that there exists apolynomial g ( y ) ∈ R [ y ] \ R dividing P ( y ) in R [ y ], a contradiction with the assumption that P ( y )is irreducible over R . Therefore, P ( x λ ) is irreducible over K ( x , . . . , x n − ).6or the second equivalence, by Gauß’ lemma, one easily sees that P ( x λ ) is irreducible over K ( x , . . . , x n − ) if and only if prim x n ( P ( x λ )) is irreducible over R [ x , . . . , x n − ]. It thus amountsto showing the equivalence between the irreducibility of prim x n ( P ( x λ )) over R [ x , . . . , x n − ] andits irreducibility over R . The direction from R to R [ x , . . . , x n − ] is trivial. In order to see theconverse, notice that any nonconstant factor of prim x n ( P ( x λ )) can only belong to R [ x , . . . , x n − ]since prim x n ( P ( x λ )) is irreducible over R [ x , . . . , x n − ]. On the other hand, the existence of anysuch a nonconstant factor would contradict with the primitivity of prim x n ( P ( x λ )) with respect to x n . Accordingly, prim x n ( P ( x λ )) must be irreducible over R . Lemma 2.3.
Let p ∈ R [ x , . . . , x n ] and λ ∈ Z n with gcd( λ , . . . , λ n ) = , λ , . . . , λ n − not allzero and λ n > . Then λ is a q-integer linear type of p if and only if there exists a q-primitivepolynomial P ( y ) ∈ K [ y ] such that P ( x λ ) divides p in K ( x , . . . , x n − )[ x n ] .Proof. The necessity is readily seen from the definition of q -integer linear types. In order to showthe su ffi ciency, it amounts to proving that there exists an irreducible polynomial in R [ x , . . . , x n ]of the form x α f ( x λ ) for α ∈ N n and f ( y ) ∈ R [ y ] dividing p in R [ x , . . . , x n ]. Let f ( y ) ∈ R [ y ] bea primitive irreducible factor of P ( y ). Since P ( y ) is q -primitive, we know that f ( y ) is q -primitiveas well. Hence prim x n ( f ( x λ )) = x α f ( x λ ) for some α ∈ N n with α n =
0. By Lemma 2.2 (ii),prim x n ( f ( x λ )) is irreducible over R . Because P ( x λ ) divides p in K ( x , . . . , x n − )[ x n ], so does f ( x λ ). It then follows that prim x n ( f ( x λ )) divides p in R [ x , . . . , x n ]. In other words, x α f ( x λ ) is anirreducible factor of p in R [ x , . . . , x n ]. Therefore, λ is a q -integer linear type of p by definition.This concludes the proof.
3. The bivariate case
Before turning to the general multivariate case, we first consider the simpler yet importantsubcase of bivariate polynomials and present a fast algorithm for computing the q -integer lineardecomposition in this context. As we will see shortly, this subcase serves as a nice illustration ofthe main idea of the first approach developed in the next section. Moreover, this subcase can beused as a base case, yielding a second, iterative approach for handling the general multivariatecase. This is discussed in Section 5. In this section, we write ( x , y ) for the two variables ( x , x ).Let p be a polynomial in R [ x , y ]. As univariate polynomials are q -integer linear, we mayassume without loss of generality that the given polynomial p is nonconstant and it is primitivewith respect to x as well as with respect to y , or equivalently, cont x ( p ) = cont y ( p ) =
1. With thisset-up, p admits the q -integer linear decomposition of the particular form p = x α P m Y i = P i ( x λ i y µ i ) , (3.1)where α, m , µ i ∈ N , λ i ∈ Z , P ∈ R [ x , y ] and P i ( y ) ∈ R [ y ] with(1) P being primitive and having only non- q -integer linear factors except for constants;(2) the P i ( y ) being nonconstant and q -primitive;(3) the ( λ i , µ i ) being distinct q -integer linear types and all λ i µ i , q -integer linear types of p . One way to proceed, according to (Le, 2001), is to compute the resultant of p and its q -shift p ( qx , q r y ) with respect to x for an indeterminate r , and then obtain the candidates ( λ, µ ) by findingpossible rationals r = − λ/µ that make the resultant zero. We show below how the problem canbe reduced to an interesting and easy geometry task, avoiding computing any resultant.Observe that, unlike the ordinary shift case, all q -integer linear types ( λ i , µ i ) in (3.1) appearas the exponent vectors, instead of as the coe ffi cients, of p . This observation implies that thehomogeneous polynomial technique from (Giesbrecht et al., 2019a) will not work in the q -case.On the other hand, it also suggests that supports of polynomials may play a role similar to ho-mogeneous polynomials in the ordinary shift case. This is exactly the key idea we are describingnow. Lemma 3.1.
Let p ∈ R [ x , y ] \ R with cont x ( p ) = cont y ( p ) = and admitting (3.1) . Then for any ℓ ∈ N with ≤ ℓ ≤ m and for any pair ( i , j ) ∈ supp( p ) , there exists another pair (˜ i , ˜ j ) ∈ supp( p ) such that ( i , j ) = (˜ i + λ ℓ k , ˜ j + µ ℓ k ) for k ∈ Z \ { } , that is, ( i − ˜ i ) / ( j − ˜ j ) = λ ℓ /µ ℓ .Proof. There is nothing to show when m =
0. Assume that m >
0. It su ffi ces to show theassertion for ℓ = m and then the lemma follows by the symmetry.Let p ∗ = x α P Q m − i = P i ( x λ i y µ i ). Then p ∗ ∈ R [ x , y ] \ { } as p ,
0. By (3.1), p = p ∗ · P m ( x λ m y µ m ) . (3.2)By assumption, supp( p ) , {} . Let ( i , j ) ∈ supp( p ). It follows from (3.2) that there is ( i ∗ , j ∗ ) ∈ supp( p ∗ ) and k ∗ ∈ supp( P m ( y )) such that ( i , j ) = ( i ∗ + λ m k ∗ , j ∗ + µ m k ∗ ). Now consider the set S ( i ∗ , j ∗ ) = { ( i ′ , j ′ ) ∈ supp( p ∗ ) | ( i ′ , j ′ ) = ( i ∗ + λ m k ′ , j ∗ + µ m k ′ ) for some k ′ ∈ Z } . Then there exist p ∗ , p ∗ ∈ R [ x , y ] with supp( p ∗ ) = S ( i ∗ , j ∗ ) and supp( p ∗ ) = supp( p ∗ ) \ S ( i ∗ , j ∗ ) suchthat p ∗ = p ∗ + p ∗ . It is evident that ( i ∗ , j ∗ ) ∈ S ( i ∗ , j ∗ ) . Thus S ( i ∗ , j ∗ ) is a nonempty set in N and then p ∗ is nonzero. Let ( α ∗ , β ∗ ) be the minimal element in S ( i ∗ , j ∗ ) . One then sees from the definitionof S ( i ∗ , j ∗ ) that any its element can be written as ( α ∗ + λ m k ′ , β ∗ + µ m k ′ ) for k ′ ∈ N , or equivalently,every monomial appearing in p ∗ has the form x α ∗ + λ m k ′ y β ∗ + µ m k ′ for k ′ ∈ N . It then follows thatthere exists a nonzero univariate polynomial P ∗ ( y ) ∈ R [ y ] such that p ∗ = x α ∗ y β ∗ P ∗ ( x λ m y µ m ).On the other hand, by noticing that for any ( i ′ , j ′ ) ∈ supp( p ∗ ) = supp( p ∗ ) \ S ( i ∗ , j ∗ ) , we have( i ′ , j ′ ) , ( i ∗ + λ m k ′ , j ∗ + µ m k ′ ) for all k ′ ∈ Z . Hence, p can be decomposed as p = f + g ,where f = p ∗ P m ( x λ m y µ m ) and g = p ∗ P m ( x λ m y µ m ) with supp( f ) ∩ supp( g ) = {} . As a consequence,supp( p ) = supp( f ) ⊎ supp( g ). Since ( i ∗ , j ∗ ) ∈ S ( i ∗ , j ∗ ) , we have ( i , j ) ∈ supp( f ). Notice that p ∗ = x α ∗ y β ∗ P ∗ ( x λ m y µ m ). So f = x α ∗ y β ∗ ˜ P ( x λ m y µ m ) with ˜ P ( y ) = P ∗ ( y ) P m ( y ) ∈ R [ y ] \ { } . Then thereexists k ∈ supp( ˜ P ( y )) such that ( i , j ) = ( α ∗ + λ m k , β ∗ + µ m k ). Since P m ( y ) is nonconstant and q -primitive, it has more than one monomial, then so does ˜ P ( y ). This implies that there is anotherelement ˜ k ∈ supp( ˜ P ( y )) with ˜ k , k . Let (˜ i , ˜ j ) = ( α ∗ + λ m ˜ k , β ∗ + µ m ˜ k ). Then the assertion for ℓ = m follows by the observation that (˜ i , ˜ j ) ∈ supp( f ) ⊂ supp( p ).The above lemma suggests a simple geometric way to find candidates for all q -integer lineartypes of the given polynomial p . For each monomial x i y j appearing in p plot a point at thecoordinate ( j , i ) on the ( y , x )-plane with y -axis the horizontal axis and x -axis the vertical axis, andthen determine the lower convex hull C of all these points (namely the convex hull of all verticalrays starting from these points and continuing upwards). We call C the Newton polygon of p .8ne can refer to (Walker, 1950, Chaper IV) for more information about the Newton polygon.The following basic properties are then geometrically evident. Lemma 3.2.
With the assumptions of Lemma 3.1, let C be the Newton polygon of p. Then(i) for each pair ( i , j ) ∈ supp( p ) , the point ( j , i ) is contained in C.(ii) for every vertex ( j , i ) of C, the pair ( i , j ) belongs to supp( p ) .(iii) with all points ( j , i ) with ( i , j ) ∈ supp( p ) plotted on the plane, the leftmost (resp. rightmost)vertex of C is one of the leftmost (resp. rightmost) points.(iv) slopes of the non-vertical edges of C increase from left to right.(v) each non-vertical edge of C has the lowest possible slope, compared with those line seg-ments inside C which start at the left vertex of the edge and end at any point on the right. The key feature of the Newton polygon of p is that the slopes of its non-vertical edges provideall possible choices for the rational numbers λ i /µ i , and thus give all candidates for the q -integerlinear types ( λ i , µ i ). Observe that all λ i µ i ,
0. Hence we consider only non-horizontal andnon-vertical edges of C , namely those which have nonzero and finite slopes. Proposition 3.3.
With the assumptions of Lemma 3.1, let C be the Newton polygon of p. Thenthe slopes of non-horizontal and non-vertical edges of C constitute a superset of all q-integerlinear types of p. Moreover, with s ∈ N denoting the cardinality of supp( p ) , this superset has nomore than s − elements in total.Proof. There is nothing to show when m = m >
0. It su ffi ces to showthat there exists one edge of C whose slope is equal to λ m /µ m , with the rest following from thesymmetry and the observation that all λ i /µ i are nonzero and finite.Suppose that none of the edges of C has the slope λ m /µ m . By assumption, one sees that p hasat least two monomials of distinct powers in y . Thus C has at least one non-vertical edge. Let A , . . . , A t + with t ∈ N \ { } be all the vertices of C from left to right, and for each 1 ≤ ℓ ≤ t let α ℓ denote the slope of the edge connecting A ℓ and A ℓ + (as visualized in Figure 1). Then α < · · · < α t < ∞ by Lemma 3.2 (iv) and we have the following case distinction. Case 1. λ m /µ m < α . Let ( j , i ) be the coordinate of the leftmost vertex A . Then ( i , j ) ∈ supp( p )by Lemma 3.2 (ii). It follows from Lemma 3.1 that there is (˜ i , ˜ j ) ∈ supp( p ) such that ˜ j , j and( i − ˜ i ) / ( j − ˜ j ) = λ m /µ m < α . Notice that C is convex with A the leftmost vertex. By Lemma 3.2(iii), the point ( ˜ j , ˜ i ) lies outside of C , a contradiction with Lemma 3.2 (i). Case 2. α ℓ < λ m /µ m < α ℓ + for some ℓ with 1 ≤ ℓ < t . Let ( j , i ) be the coordinate of the vertex A ℓ . A similar argument as Case 1 leads to the contradiction that there exists a point ( ˜ j , ˜ i ) with(˜ i , ˜ j ) ∈ supp( p ) lying outside of C . Case 3. α t < λ m /µ m . Let ( j , i ) be the coordinate of the rightmost vertex A t + . A similar argumentas Case 1 leads to the contradiction that there exists a point ( ˜ j , ˜ i ) with (˜ i , ˜ j ) ∈ supp( p ) lyingoutside of C .The above discussions provide the following necessary condition for p to be a q -integer linearpolynomial. 9 m /µ m λ m /µ m λ m /µ m xO y • • • • • • • • • • • • •• • • • • • •• A α A A ℓ α ℓ A ℓ + α ℓ + A ℓ + A t α t A t + C Figure 1: An illustration of the Newton polygon C of p in Proposition 3.3. Proposition 3.4.
With the assumptions of Lemma 3.1, further assume that p is q-integer linear.Then there is only one monomial appearing in p of lowest or highest degree in y; a symmetry ar-gument holds in terms of x. Geometrically speaking, this is equivalent to say that, with all points ( j , i ) corresponding to monomials x i y j appearing in p plotted on the ( y , x ) -plane, there is onlyone leftmost point, one rightmost point, one lowest point and one highest point. Consequently,the Newton polygon of p does not have horizontal edges.Proof. We only show the assertion for monomials of lowest degree in y , corresponding to left-most points on the ( y , x )-plane, with the other three assertions following by similar arguments.In this respect, it is su ffi cient to show that if p contains two distinct monomials of lowestdegree in y then so does its non- q -integer linear part, namely P in (3.1), because then we knowthat P is not a unit in R , and thus obtain the contradiction that p is not q -integer linear.We proceed by induction on the number m in (3.1). The assertion is evident when m = p = x α P by (3.1). Assume that m > m −
1, namely for anypolynomial in R [ x , y ] of form (3.1) with m replaced by m −
1, if it has two di ff erent monomialsof lowest degree in y then so does its non- q -integer linear part. Let j ∈ N be the lowest degreeof p in y and suppose that x i y j and x i y j with i , i in N are two monomials appearing in p . Let p ∗ = x α P Q m − i = P i ( x λ i y µ i ) so that (3.2) holds. It amounts to showing that x i y j and x i y j both appear in p ∗ , which, along with a subsequent application of the induction hypothesisto p ∗ , concludes the proof. We see from (3.2) that there exist ( i ∗ , j ∗ ) , ( i ∗ , j ∗ ) ∈ supp( p ∗ ) and k , k ∈ supp( P m ( y )) such that j = j ∗ + µ m k = j ∗ + µ m k , i = i ∗ + λ m k , i = i ∗ + λ m k . Since µ m > P m ( y ) is q -primitive, we derive from the minimality of j that k = k = j ∗ = j ∗ = j . Therefore, i ∗ = i and i ∗ = i . This implies that ( i , j ) , ( i , j ) ∈ supp( p ∗ ),that is, x i y j and x i y j indeed both appear in p ∗ .With candidates for the q -integer linear types ( λ i , µ i ) at hand, we are now able to find thecorresponding univariate polynomials P i ( y ) based on a q -counterpart of (Giesbrecht et al., 2019a,Proposition 3.2). 10 roposition 3.5. With the assumptions of Lemma 3.1, let ( λ, µ ) ∈ Z with gcd( λ, µ ) = , λ , and µ > . Let P ∗ ( y ) ∈ R [ y ] be the content of the numerator of p ( x µ , yx − λ ) with respect to x. Then ( λ, µ ) is a q-integer linear type of p if and only if P ∗ ( y ) < R . Moreover, if ( λ, µ ) is a q-integerlinear type of p then P ∗ ( y /µ ) ∈ R [ y ] and it is the corresponding univariate polynomial of thetype ( λ, µ ) .Proof. For the necessity of the first assertion, assume that ( λ, µ ) is a q -integer linear type of p .Then by (3.1), there exists an integer i with 1 ≤ i ≤ m such that ( λ i , µ i ) = ( λ, µ ), and thus p ( x µ , yx − λ ) = x µα P ( x µ , yx − λ ) m Y j = , j , i P j ( x λ j µ − λµ j y µ j ) · P i ( y µ ) ∈ R [ x , x − , y ] . Notice that µ > P i ( y ) ∈ R [ y ] \ R . Therefore one sees from the above equation that P i ( y µ )divides all coe ffi cients of p ( x µ , yx − λ ) with respect to x , giving rise to a desired nontrivial commondivisor of P ∗ ( y ) in R [ y ].To show the su ffi ciency, let f ( y ) ∈ K [ y ] \ K be a monic irreducible factor of P ∗ ( y ). Then f ( y )divides p ( x µ , yx − λ ) in K ( x )[ y ]. Substituting y by yx λ and x by x /µ yields that f ( yx λ/µ ) divides p in K ( x )[ y ] with K ( x ) the algebraic closure of the field K ( x ). This implies that f ( y ) , y , for,otherwise, we would have y divides p in K ( x )[ y ] and then p ( x , =
0, a contradiction with theprimitivity of p with respect to y . Observe that f ( y ) is a monic irreducible polynomial in K [ y ] ofpositive degree. Let r ∈ K be a root of f ( y ). Then r , f ( y ) is its minimal polynomial in K [ y ]. Meanwhile, it follows from the divisibility of p by f ( yx λ/µ ) in K ( x )[ y ] that y = rx − λ/µ is aroot of p in K ( x ), that is, p ( x , rx − λ/µ ) =
0. Let P ( y ) ∈ K [ y ] be the minimal polynomial of r µ . Thendeg( P ( y )) > P (0) , r ,
0. Also, by Lemma 2.2 (ii), P ( x λ y µ ) is irreducible over K ( x ).A simple calculation verifies that P ( x λ y µ ) also vanishes at y = rx − λ/µ . Therefore, P ( x λ y µ ), uponmaking it monic with respect to y , gives rise to the minimal polynomial of y = rx − λ/µ in K ( x )[ y ].One then derives from p ( x , rx − λ/µ ) = P ( x λ y µ ) divides p in K ( x )[ y ]. By Lemma 2.3, ( λ, µ )is a q -integer linear type of p .From the preceding paragraph, one can additionally conclude that f ( y ) divides P ( y µ ) in K [ y ],and then by Gauß’ lemma, prim y ( f ( y )) divides prim y ( P ( y µ )) in R [ y ]. Indeed, by the definition of P ( y ), one sees that P ( r µ ) =
0. This in turn means that r is a root of P ( y µ ) in K . Since P ( y µ ) ∈ K [ y ]and f ( y ) is the minimal polynomial of r in K [ y ], one confirms that f ( y ) divides P ( y µ ) in K [ y ].It remains to show the last assertion. Notice that p is primitive with respect to y , in particular p ( x , ,
0. Then P ∗ ( y ) is q -primitive. Since ( λ, µ ) is a q -integer linear type of p , then ( λ, µ ) = ( λ i , µ i ) for some integer i with 1 ≤ i ≤ m . Thus it amounts to showing that P ∗ ( y ) = P i ( y µ ),which in turn su ffi ces to prove that P i ( y µ ) divides P ∗ ( y ) and conversely in R [ y ], because then theassertion directly follows from the fact that both P ∗ ( y ) and P i ( y ) are q -primitive.As in the proof of the necessity, one readily sees from (3.1) that P i ( y µ ) divides P ∗ ( y ) in R [ y ]. In order to show the converse, by noticing that P ∗ ( y ) is nonconstant as deg( P i ( y )) >
0, let f ( y ) ∈ R [ y ] be a primitive irreducible factor of P ∗ ( y ) of positive degree in y . From the proofof the su ffi ciency, one concludes that there exists P ( y ) ∈ K [ y ] with P (0) , P ( x λ y µ )divides p in K ( x )[ y ] and f ( y ) divides prim y ( P ( y µ )) in R [ y ]. By (3.1), P ( x λ y µ ) divides P i ( x λ y µ )in K ( x )[ y ]. According to Lemma 2.2 (i), there exists c ∈ K , β ∈ Z and a factor g ( y ) ∈ R [ y ] of P i ( y ) such that P ( x λ y µ ) = cx β g ( x λ y µ ). Since P (0) ,
0, we know that β =
0, that is, P ( y ) = cg ( y ).This implies that prim y ( P ( y )) divides P i ( y ) in R [ y ]. Therefore, f ( y ) divides P i ( y µ ) in R [ y ]. By thearbitrariness of f ( y ), one obtains that P ∗ ( y ) divides P i ( y µ ) in R [ y ]. This completes the proof.11ote that the previous proposition also guarantees that any fake candidate for q -integer lineartypes can be easily recognized by a content computation. Also note that fake candidates onlycome from the polynomial P in (3.1). Therefore, any occurrence of a trivial content immediatelyindicates that the given polynomial is not q -integer linear.Putting things together, we obtain a new algorithm for computing the q -integer linear decom-position of a bivariate polynomial. For compatibility with later algorithms, we take a nonconstantand primitive polynomial as input. BivariateQILD.
Given a nonconstant and primitive polynomial p ∈ R [ x , y ], compute its q -integer linear decomposition.1. Set f ( y ) = cont x ( p ), P = prim x ( p ), m =
0, and set β to be the lowest degree of f ( y )with respect to y . Update f ( y ) = f ( y ) / y β . If f ( y ) , m = m +
1, and set( λ m , µ m ) = (0 ,
1) and P m ( y ) = f ( y ).2. Set f ( x ) = cont y ( P ) and set α to be the lowest degree of f ( x ) with respect to x . Update P = prim y ( P ) and f ( x ) = f ( x ) / x α . If f ( x ) , m = m +
1, and set( λ m , µ m ) = (1 ,
0) and P m ( y ) = f ( y ).3. If deg( P ) = cx α y β P Q mi = P i ( x λ i y µ i ).4. Find the support supp( P ) = { ( i , j ) , . . . , ( i s , j s ) } ⊂ N with ( i , j ) < · · · < ( i s , j s ).5. Set Λ = {} and k = j k , j s do5.1 Find the maximum ℓ ∈ { k + , . . . , s } such that j ℓ , j k and r = i ℓ − i k j ℓ − j k is minimal.5.2 Update k = ℓ .5.3 If r , Λ = Λ ∪ { ( λ, µ ) ∈ Z | λ/µ = r , gcd( λ, µ ) = µ > } .6. For ( λ, µ ) in Λ do6.1 Set P ∗ ( y ) to be the content of the numerator of P ( x µ , yx − λ ) with respect to x .6.2 If deg( P ∗ ( y )) > m = m +
1, ( λ m , µ m ) = ( λ, µ ), P m ( y ) = P ∗ ( y /µ ).Set f , g ∈ R [ x , y ] to be the numerator and denominator of P m ( x λ y µ ), andupdate P = P / f and α = α + deg x ( g ).7. Return cx α y β P Q mi = P i ( x λ i y µ i ). Theorem 3.6.
Let p ∈ R [ x , y ] be a valid input of the algorithm BivariateQILD . Then the algo-rithm terminates and correctly computes the q-integer linear decomposition of p.Proof.
Notice that the integer k increases every time the algorithm passes through Steps 5.1-5.3.Thus there are only finitely many iterations happening in Step 5, and then the cardinality of theset Λ is finite. Therefore, the algorithm terminates.In order to prove the correctness, it su ffi ces to show that the set Λ constructed in Step 6 isin one-to-one correspondence with the set comprised of all slopes of non-horizontal and non-vertical edges of the Newton polygon of p , with the rest following from Propositions 3.3 and 3.5.This in turn is evident by Lemma 3.2 (v). 12 emark 3.7. Step 5 of the above algorithm presents a naive way to construct the Newton polygonof a bivariate polynomial. Fast algorithms are available (cf. (Graham, 1972; Anderson, 1978)).Nevertheless, the cost of this step will not dominate the total cost. Rather than making e ff ort togain a slight improvement, we prefer to present the computation in such a way that it extends toa higher dimensional space in a straightforward manner (see Section 4). Remark 3.8.
If one is merely interested in the q-integer linearity of the input polynomial p ∈ R [ x , y ] , rather than the full q-integer linear decomposition, then our algorithm can be easilymodified: any of the following conditions will trigger the adapted algorithm to terminate early,returning that p is not q-integer linear. • (Proposition 3.4) In Step 4, the support supp( p ) has more than one element whose eitherentry has the extremum value. Note that this includes the case when r = in Step 5.3. • (Proposition 3.5) In Step 6.2, the case of deg( P ∗ ( y )) = happens, that is, the candidate ( λ, µ ) currently under investigation is fake. • (Definition 2.1) In Step 7, we have deg( P ) > .3.1. Complexity analysis of the bivariate algorithm Although our algorithms work in more general UFDs, we confine our complexity analysisto the case of integer (Laurent) polynomials, that is, when D is the ring of integers Z and then R is equal to Z [ q , q − ]. Here q can be viewed as a variable in addition to x , . . . , x n . Note thatoperations in Z [ q , q − ] can be easily transferred to those in Z [ q ] with a negligible cost. The costis given in terms of number of word operations used so that growth of coe ffi cients comes intoplay. Recall that the word length of a nonzero integer a ∈ Z is defined as O(log | a | ). In thispaper, all complexity is analyzed in terms of O-estimates (or O ∼ -estimates) for classical or fastarithmetic, where the soft-Oh notation “O ∼ ” is basically “O” but suppressing logarithmic factors(see (von zur Gathen and Gerhard, 2013, Definition 25.8) for a precise definition).Throughout this paper, we define the max-norm || p || ∞ of a Laurent polynomial p ∈ Z [ q , q − ]as the maximum absolute value of its coe ffi cients with respect to q , and the max-norm || p || ∞ of a polynomial p = P i ∈ N p i ,..., i n x i ∈ Z [ q , q − ][ x , . . . , x n ] as max i ∈ N {|| p i ,..., i n || ∞ } . The GCDcomputation is fundamental for our algorithms. Before analyzing the algorithm, let us recallsome useful complexity results on GCD computation. Lemma 3.9 ((Gel ′ fond, 1960, Page 135-139)) . Let p , . . . , p m ∈ Z [ x , . . . , x n ] . Let p = p · · · p m and let d i = deg x i ( p ) for all i = , . . . , n. Then || p || ∞ · · · || p m || ∞ ≤ e d + ··· + d n || p || ∞ , where e is the Euler constant. Therefore, log || p i || ∞ ∈ O( d + log || p || ∞ ) for all i = , . . . , n withd = max { d , . . . , d n } . Note that when n = d , which, however, leads to the same order of magnitude for word lengths of the max-norms. Lemma 3.10.
Let f , g ∈ Z [ x , . . . , x n ] with deg x i ( f ) , deg x i ( g ) ≤ d i for all i = , . . . , n and || f || ∞ , || g || ∞ ≤ β . Let d = max { d , . . . , d n } and D n = d · · · d n . Then computing gcd( f , g ) over Z takes O ∼ ( d D n + D n log β ) word operations with classical arithmetic and O ∼ ( dD n + D n log β ) with fast arithmetic. roof. We proceed to compute gcd( f , g ) by a small prime modular algorithm. By Lemma 3.9, || gcd( f , g ) || ∞ ≤ e d + ··· + d n β ≤ e nd β = B with e the Euler constant. Then log B ∈ O( nd + log β ).It is su ffi cient to choose ⌈ log (2 B + ⌉ primes, each of size O(log B ). For every chosen prime h , we then reduce all coe ffi cients of f and g modulo h , using O( D n log β log h ) word operationswith classical arithmetic, and compute gcd( f h , g h ) with f h = f mod h and g h = g mod h . Thedesired gcd( f , g ) can be recovered by a final application of the Chinese remainder theorem, whichtakes O ∼ ( d D n + D n log β ) word operations with classical arithmetic and O ∼ ( dD n + D n log β )with fast arithmetic as each log h ∈ O(log log B ) and n ∈ O(log D n ). Now it remains to count thenumber of arithmetic operations, denoted by G h ( n , d , D n ), used by the gcd computation in thefield Z h for each prime h , with the rest following by the fact that each operation of these takesO(log h ) word operations with classical arithmetic and O ∼ (log h ) with fast arithmetic.For each prime h , we compute gcd( f h , g h ) with f h = f mod h and g h = g mod h by anevaluation-interpolation scheme: evaluate coe ffi cients of f h , g h with respect to x , . . . , x n − at d n points from Z h for x n ; compute d n GCDs over Z h of two ( n − d , . . . , d n − in x , . . . , x n − , respectively; recover the final GCD by interpolation. Noticethat there are at most d · · · d n − monomials in x , . . . , x n − appearing in each of f h and g h . Theevaluation and interpolation steps take O( d n D n ) arithmetic operations in Z h with classical arith-metic and O ∼ ( D n ) with fast arithmetic. The second step uses O( d n G h ( n − , d , D n − )) arithmeticoperations in Z h , where D n − = d · · · d n − . Thus we obtain the recurrence relationO( G h ( n , d , D n )) ⊂ O( d n D n ) + O( d n G h ( n − , d , D n − ))with classical arithmetic andO( G h ( n , d , D n )) ⊂ O ∼ ( D n ) + O( d n G h ( n − , d , D n − ))with fast arithmetic. From the initial condition that G h (1 , d , d ) is in O( d ) with classical arith-metic and in O ∼ ( d ) with fast arithmetic, one concludes that G h ( n , d , D n ) is in O( ndD n ) withclassical arithmetic and in O ∼ ( D n ) with fast arithmetic.We are now ready to present the complexity of our bivariate algorithm. In order to make itready to use in the analysis of our second approach in Section 5, we analyze the cost in the caseof R = Z [ q , q − , z , . . . , z v ], where v ∈ N is arbitrary but fixed and the z i are additional parametersindependent of q , x , y . Theorem 3.11.
Let p ∈ Z [ q , q − , z , . . . , z v ][ x , y ] . Assume that the numerator and denominatorof p have maximum degree d in each variable from { q , z , . . . , z v , x , y } separately and let || p || ∞ = β . Then the algorithm BivariateQILD computes the q-integer linear decomposition of p over Z [ q , q − , z , . . . , z v ] using O ∼ ( d v + + d v + log β ) word operations with classical arithmetic and O ∼ ( d v + + d v + log β ) with fast arithmetic. In particular, when v = the algorithm uses O ∼ ( d + d log β ) word operations with classical arithmetic and O ∼ ( d + d log β ) with fast arithmetic.Proof. In Step 1, finding the content f ( y ) amounts to computing a GCD of at most ( d + Z [ q , z , . . . , z v , y ] of maximum degrees at most d in each variable separately andmax-norms at most β . Thus by Lemma 3.10, this step takes O ∼ ( d v + + d v + log β ) word operationswith classical arithmetic and O ∼ ( d v + + d v + log β ) with fast arithmetic. The same cost applies toStep 2. Step 3 is the trivial case and takes no word operations. Step 4 takes linear time in thecardinality s of the support, namely O( d ) word operations. In Step 5, it is readily seen that eachiteration of the loop takes O(( s − k ) log( s − k ) log d ) word operations with classical arithmetic,14ielding the total cost O( d log d ) word operations with classical arithmetic as there are at most s ≤ ( d + iterations. In Step 6, for each ( λ, µ ) ∈ Λ , a direct calculation shows that P ( x µ , yx − λ )and P have the same degree in y which is at most d , the same max-norm which is of word lengthO(( v + d + log β ) and the same number of nonzero monomials in x , y appearing which is atmost ( d + . Thus by Lemma 3.10, Step 6.1 takes O ∼ ( d v + + d v + log β ) word operations withclassical arithmetic and O ∼ ( d v + + d v + log β ) with fast arithmetic, which dominates the cost forStep 6.2. Note that we can easily expand P ( x µ , yx − λ ) within the allowed costs. Since there areat most s − ≤ ( d + − Λ , the claimed cost then follows. Remark 3.12.
If one finds a multivariate version of the algorithm of Conflitti (2003), then thecomplexity can be further improved.
4. The first approach for the multivariate case
In this and the next section, we will deal with multivariate polynomials having more than twovariables. We propose two approaches to compute the q -integer linear decomposition of such apolynomial. The first one extends the bivariate algorithm introduced in the preceding sectionto the multivariate case in a straightforward manner, which is discussed in the present section.The second approach takes the bivariate algorithm as a base case and proceeds in an iterativefashion by tackling two variables at each iteration. This is explored in the next section. Unlikethe ordinary shift case in (Giesbrecht et al., 2019a), both approaches appear to perform well inpractice as we shall see in Section 7.Let us now start with the first approach. This method follows exactly the same pattern asthe bivariate algorithm given in Section 3: first find all possible candidates for q -integer lineartypes of the given polynomial by means of geometry; then extract the corresponding univariatepolynomials for the types one by one via content computations. All results presented below canbe shown along almost the same lines as in the bivariate case but from a multivariate point ofview, so we will omit the proofs.Let p ∈ R [ x , . . . , x n ] and assume that it admits the q -integer linear decomposition (2.1).By iteratively removing the content of p with respect to each variable from { x , . . . , x n } andrecursively computing the q -integer linear decomposition of every removed content which is an( n − R , we may further assume without loss of generality that p isnonconstant and primitive with respect to each variable from { x , . . . , x n } . This implies that inthe equation (2.1), we have c = α n = λ i has zero entries.In order to find the candidates, one shows the following multivariate version of Lemma 3.1. Lemma 4.1.
Let p ∈ R [ x , . . . , x n ] \ R with cont x ( p ) = · · · = cont x n ( p ) = . Then for any j ∈ N with ≤ j ≤ m and for any i ∈ supp( p ) , there exists another n-vector ˜ i ∈ supp( p ) such that i = ˜ i + k · λ j for some k ∈ Z \ { } . By the same reasoning as in the bivariate case, the problem is reduced to computing “slopes”of edges of the
Newton polytope associated to the given polynomial p . By the Newton polytopeassociated to p , we mean the lower convex hull (in view of the x -axis) in the ( x , . . . , x n )-spaceof all points corresponding to vectors in supp( p ), whose edges are comprised of line segmentsconnecting the minimal vector to the maximal one in supp( p ) in terms of the total ordering “ < ”on N n . And the “slope” of the line connecting two points i , j ∈ N n in the ( x , . . . , x n )-space isreferred to as the rational vector ( i − j i n − j n , . . . , i n − − j n − i n − j n , ∈ Q n .15 roposition 4.2. With the assumptions of Lemma 4.1, the finite “slopes” of edges of the Newtonpolytope associated to p which have no zero entries constitute a superset of all q-integer lineartypes of p. Moreover, with s ∈ N denoting the cardinality of supp( p ) , this superset has no morethan s − elements in total. Similar to Proposition 3.4, we obtain a necessary condition for p to be q -integer linear. Proposition 4.3.
With the assumptions of Lemma 4.1, if p is q-integer linear, then the extremumvalue of every variable from { x , . . . , x n } amongst the vectors of supp( p ) can only be attainedonce. As a consequence, all “slopes” of edges of the Newton polytope associated to p have nozero entries. In order to generalize Step 5 of the bivariate algorithm to the multivariate case, we extendthe total ordering “ < ” on N n to Q n as follows. Let i , j ∈ Q n , we say that i < j if the rightmostnonzero entry of the vector di ff erence i − j ∈ Q n is negative. Due to its convexity, the constructionof the Newton polytope associated to p can be visualized in the following steps:(i) plot a point at the coordinate i in the ( x , . . . , x n )-space for every vector i ∈ supp( p );(ii) take the point corresponding to the minimal element in supp( p ) and denote it by A ;(iii) draw a line segment from A to every point strictly on its right in view of the x n -axis;(iv) find the line segment of lowest possible “slope”, yielding one edge of the Newton polytope;(v) update A to be the other end point of this edge;(vi) repeat (iii)-(v) until A reaches the rightmost point in view of the x n -axis.The above description, together with Proposition 4.2, provides a natural way to compute can-didates for q -integer linear types of p . In analogy to Proposition 3.5, we have the followingmultivariate criterion for detecting the genuineness of each candidate, along with the computa-tion of the corresponding univariate polynomial when the answer is a ffi rmative. Proposition 4.4.
With the assumptions of Lemma 4.1, let λ ∈ Z n with gcd( λ , . . . , λ n ) = , λ , . . . , λ n − not all zero and λ n > . Let P ∗ ( y ) ∈ R [ y ] be the content of the numerator ofp ( x λ n , . . . , x λ n n − , yx − λ · · · x − λ n − n − ) with respect to x , . . . , x n − . Then λ is a q-integer linear type ofp if and only if P ∗ ( y ) < R. Moreover, if λ is a q-integer linear type of p then P ∗ ( y /λ n ) ∈ R [ y ] andit is the corresponding univariate polynomial of the type λ . Assembling everything together yields our first approach.
MultivariateQILD . Given a polynomial p ∈ R [ x , . . . , x n ], compute its q -integer linear decom-position.1. If p ∈ R then set c = p ; and return c .2. Set c = cont( p ) and f = prim( p ). If supp( f ) is a singleton then set α to be the only elementand update c = c f / x α ; and return c x α .3. If n = α to be the lowest degree of f with respect to x , m = λ m = P m ( y ) = f ( y ) / y α ; and return c x α Q mi = P i ( x λ i ).16. Set α = , P = m = i = , . . . , n do4.1 Set g = cont x i ( f ), and update f = prim x i ( f ).4.2 If g , g ∈ R [ x , . . . , x i − , x i + , . . . , x n ],returning g = x ˜ α · · · x ˜ α i − i − x ˜ α i + i + · · · x ˜ α n n ˜ P m Y j = ˜ P j ( x ˜ λ j · · · x ˜ λ j , i − i − x ˜ λ j , i + i + · · · x ˜ λ jn n ) , update α = α + ( ˜ α , . . . , ˜ α i − , , ˜ α i + , . . . , ˜ α n ), P = P ˜ P , and for j = , . . . , ˜ m iteratively update m = m + λ m = ( ˜ λ j , . . . , ˜ λ j , i − , , ˜ λ j , i + , . . . , ˜ λ jn ), P m ( y ) = ˜ P j ( y ).5. If deg( f ) = c = c f ; and return c x α P Q mi = P i ( x λ i ).6. Find the support supp( f ) = { i , . . . , i s } with i < · · · < i s .7. Set Λ = {} and k = i kn , i sn do7.1 Find the maximum j ∈ { k + , . . . , s } such that i jn , i kn and r = (cid:18) i j − i k i jn − i kn , . . . , i j , n − − i k , n − i jn − i kn , (cid:19) is minimal.7.2 Update k = j .7.3 If r has no zero entries then update Λ = Λ ∪ { λ ∈ Z n | ( λ /λ n , . . . , λ n − /λ n , = r , gcd( λ , . . . , λ n ) = , λ n > } .
8. For λ in Λ do8.1 Set P ∗ ( y ) to be the content of the numerator of f ( x λ n , . . . , x λ n n − , yx − λ · · · x − λ n − n − ) withrespect to x , . . . , x n − .8.2 If deg( P ∗ ( y )) > m = m + λ m = λ , P m ( y ) = P ∗ ( y /λ n ).Set f ∗ , g ∗ ∈ R [ x , . . . , x n ] to be the numerator and denominator of P m ( x λ ),and update f = f / f ∗ and α i = α i + deg x i ( g ∗ ) for i = , . . . , n − f ) > P = P f else update c = c f .10. Return c x α P Q mi = P i ( x λ i ). Theorem 4.5.
Let p ∈ R [ x , . . . , x n ] . Then the algorithm MultivariateQILD terminates andcorrectly computes the q-integer linear decomposition of p.Proof. This is evident by Propositions 4.2, 4.4 and the discussions in between.
Remark 4.6.
If one is merely interested in the q-integer linearity of the input polynomial p ∈ R [ x , . . . , x n ] , rather than the full q-integer linear decomposition, then the above algorithm canbe easily modified. Analogous to Remark 3.8, any of the following conditions will trigger theadapted algorithm to terminate early, returning that p is not q-integer linear. • In Step 4.2, each polynomial g turns out to be non-q-integer linear. (Proposition 4.3) In Step 6, the support supp( f ) has more than one element whose eitherentry has the extremum value. Note that this includes the case when the integer vector r has a zero entry in Step 7.3. • (Proposition 4.4) In Step 8.2, the case of deg( P ∗ ( y )) = happens, that is, the candidate λ currently under investigation is fake. • (Definition 2.1) In Step 10, we have deg( P ) > . Theorem 4.7.
Let p ∈ Z [ q , q − ][ x , . . . , x n ] . Assume that the numerator and denominator ofp have maximum degree d in each variable from { q , x , . . . , x n } separately and let || p || ∞ = β .Then the algorithm MultivariateQILD computes the q-integer linear decomposition of p over Z using O ∼ ( n n d n + + n n d n + log β ) word operations with classical arithmetic and O ∼ ( n n d n + + n n d n + log β ) with fast arithmetic.Proof. Let T ( n , d , log β ) denote the number of word operations used by the algorithm applied to p . Steps 1 and 5 treat the trivial case, taking no word operations. In Step 2, finding the content c amounts to computing a GCD of at most ( d + n polynomials in Z [ q ] of degrees in q at most d and max-norms at most β . Thus by Lemma 3.10, this step takes O ∼ ( d n + + d n + log β ) wordoperations with classical arithmetic and O ∼ ( d n + + d n + log β ) with fast arithmetic. Step 3 dealswith the univariate case, yielding that the initial cost T (1 , d , log β ) is in O ∼ ( d + d log β ) withclassical arithmetic and O ∼ ( d + d log β ) with fast arithmetic.In Step 4, at each iteration of the loop, the computation of the content g and its primitive partin Step 4.1 can be done within O ∼ ( d n + + d n + log β ) word operations with classical arithmeticand O ∼ ( d n + + d n + log β ) with fast arithmetic; while Step 4.2 takes O( T ( n − , d , nd + log β )) wordoperations as g ∈ Z [ q , x , . . . , x i − , x i + , . . . , x n ] of maximum degree at most d in each variableand max-norm of word length O( nd + log β ) by Lemma 3.9. Since there are n iterations, thisstep in total takes O ∼ ( d n + + d n + log β ) word operations with classical arithmetic and O ∼ ( d n + + d n + log β ) with fast arithmetic, plus O( nT ( n − , d , nd + log β )) word operations.Step 6 takes linear time in the cardinality s of the support, namely O( d n ) word operations.In Step 7, each iteration of the loop requires O(( s − k ) log( s − k )( n −
1) log d ) word operationswith classical arithmetic, yielding the total cost O( d n n log d ) word operations with classicalarithmetic as there are at most s ≤ ( d + n iterations. In Step 8, for each λ ∈ Λ , a direct calculationshows that f ( x λ n , . . . , x λ n n − , yx − λ · · · x − λ n − n − ) and f have the same degree in y which is at most d ,the same max-norm which is of word length O( nd + log β ) and the same number of nonzeromonomials in x , . . . , x n − , y appearing which is at most ( d + n . Thus by Lemma 3.10, Step 8.1takes O ∼ ( d n + + d n + log β ) word operations with classical arithmetic and O ∼ ( d n + + d n + log β )with fast arithmetic, which dominates the cost for Step 8.2. Since there are at most s − ≤ ( d + n − Λ , this step takes O ∼ ( d n + + d n + log β ) word operations withclassical arithmetic and O ∼ ( d n + + d n + log β ) with fast arithmetic. Steps 9 and 10 both take noword operations without expanding the product.In summary, we obtain the recurrence relationO( T ( n , d , log β )) ⊂ O ∼ ( d n + + d n + log β ) + O( nT ( n − , d , nd + log β )) , along with T (1 , d , log β ) ∈ O ∼ ( d + d log β ) with classical arithmetic orO( T ( n , d , log β )) ⊂ O ∼ ( d n + + d n + log β ) + O( nT ( n − , d , nd + log β )) , along with T (1 , d , log β ) ∈ O ∼ ( d + d log β ) with fast arithmetic. The announced cost follows.18 O x • • • • • • • • • • • • • •• • • • (0 ,
13) (1 ,
12) (15 ,
0) (30 , Figure 2: Newton polygon of p in Example 4.8 when viewed as a polynomial in x , x . Example 4.8.
Consider the polynomial p ∈ Z [ q , q − ][ x , x , x , x ] of the formp = q x x x + qx x x + qx x x x + q x x x x + qx x x x + qx x x x − qx x x x − x x x x − x x x x − qx x x x − x x x x − x x x x + q x x x x + qx x x + qx x + q x x x + q x x x x + q x x x x + q x x x x + qx x x x + qx x x x − q x x x x − q x x x x − q x x x x + q x x x x + q x x x x + q x x x x (4.1) In order to compute the q-integer linear decomposition of the polynomial p over Z [ q , q − ] , thealgorithm MultivariateQILD first tries to find candidates for all possible q-integer linear typesof p. In this respect, it computes the support of p, which can be readily read out from (4.1) : n (9 , , , < (8 , , , < (8 , , , < (11 , , , < (10 , , , < (10 , , , < (5 , , , < (4 , , , < (4 , , , < (7 , , , < (6 , , , < (6 , , , < (1 , , , < (0 , , , < (0 , , , < (15 , , , < (14 , , , < (14 , , , < (3 , , , < (2 , , , < (2 , , , < (11 , , , < (10 , , , < (10 , , , < (7 , , , < (6 , , , < (6 , , , o . From above, we can already tell from Proposition 4.3 that the given polynomial p is not q-integer linear, since there are two elements of supp( p ) (namely the first two elements) attainingthe minimum value of x .By the definition of the ordering “ < ” on Q n , one readily sees that the projection of the Newtonpolytope associated to p on the ( x , x ) -plane coincides with the Newton polygon of p whenviewed as a polynomial in x , x . The latter is depicted in Figure 2, revealing that the ( x , x ) -coordinates of all its vertices are given by (0 , , (1 , , (15 , , (30 , . Notice that each of the ast three points gives rise to a unique element in supp( p ) . It then follows that the following fourpoints: (9 , , , , (8 , , , , (0 , , , , (6 , , , exhaust all vertices of the Newton polytope associated to p in the ( x , x , x , x ) -space. This thusyields three candidates for q-integer linear types of p, namely ( − , , − , , (2 , − , , , ( − , , , . A subsequent content computation for each candidate leads to the following q-integer lineardecomposition p = x x x · P · P ( x x − x x ) · P ( x − x x − x ) , (4.2) where P = qx x + x x + x x , P ( y ) = q y + qy + and P ( y ) = qy − y + q. Notethat the candidate ( − , , − , is fake, implying, once again, that the given polynomial p is notq-integer linear.
5. The second approach for the multivariate case
This section presents the second approach for computing the q -integer linear decompositionof a polynomial in an arbitrary number of variables. In order to describe it, we need a q -analogueof (Abramov and Petkovˇsek, 2002, Proposition 7). To this end, we require two technical lemmas.The first one corresponds to (Abramov and Petkovˇsek, 2002, Lemma 2) but restricted to the caseof Laurent polynomials. Lemma 5.1.
Let p ∈ R [ x , x − ] be a nonzero Laurent polynomial. If there exists a nonzerointeger a and a nonzero element c ∈ R such that p ( q a x ) = cp ( x ) , then c = q am for some m ∈ Z and p ( x ) / x m ∈ R .Proof. Notice that p is nonzero. Write p ( x ) = x m f ( x ) for m ∈ Z and f ∈ R [ x ] with f (0) ,
0. Itfollows from p ( q a x ) = cp ( x ) that q am f ( q a x ) = c f ( x ). Letting x = c = q am as f (0) ,
0. Hence f ( q a x ) = f ( x ). Let x be a nonzero element in R . Then f ( x ) ∈ R . Byinduction on k , we see that f ( q ka x ) = f ( x ) for all k ∈ N . Let g ( x ) = f ( x ) − f ( x ). Then g ∈ R [ x ]and it vanishes on { q ka x | k ∈ N } , which is an infinity set in R since the characteristic of R is zeroand q ∈ R is an indeterminate. Therefore g is the zero polynomial and thus f ( x ) = f ( x ) ∈ R .The lemma follows.Evidently, in the above lemma, the ring R can be replaced by any of its ring extensions whichis independent of the variable x . The next lemma plays the role of (Abramov and Petkovˇsek,2001, Lemma 3), or identically, (Hou, 2004, Lemma 3.3), in the q -shift setting, which describesa nice structure of q -shift invariant bivariate polynomials. Lemma 5.2.
Let p ∈ R [ x , y ] . If there exists c ∈ R and a , b ∈ Z , not both zero, such thatp ( q a x , q b y ) = cp ( x , y ) , then there is a univariate polynomial P ( y ) ∈ R [ y ] and four integers α, β, λ, µ with λ, µ not both zero such that p = x α y β P ( x λ y µ ) ; and conversely.Proof. There is nothing to show if p =
0. Now assume that p is nonzero and we adapt the ideafrom the proof of (Hou, 2004, Lemma 3.3) into the current setting. Let d = gcd( a , b ). Then d , a , b are not both zero. Thus letting λ = − b / d and µ = a / d gives us two coprimeintegers λ, µ . By B´ezout’s relation, there exist s , t ∈ Z such that s λ + t µ =
1. Now define20 ( x , y ) = p ( x µ y s , x − λ y t ). Then f is a nonzero Laurent polynomial in the ring R [ x , x − , y , y − ] and p ( x , y ) = f ( x t y − s , x λ y µ ) by B´ezout’s relation. Using p ( q a x , q b y ) = cp ( x , y ), a direct calculationshows that f ( q d x , y ) = c f ( x , y ). Since d ,
0, applying Lemma 5.1 to f with R replaced by R [ y , y − ] and extracting a nonpositive power of y yields that c = q dm for m ∈ Z and f ( x , y ) = x m y − k P ( y ) for k ∈ N and P ( y ) ∈ R [ y ]. Therefore, p = x tm − λ k y − sm − µ k P ( x λ y µ ). The assertionfollows by letting α = tm − λ k and β = − sm − µ k . The converse argument is evident by setting a = µ and b = − λ . This completes the proof.From the above lemma, we are then able to establish the fact that the problem of multivari-ate q -integer linearity is actually made up of a collection of subproblems of bivariate q -integerlinearity. Proposition 5.3.
Let p ∈ R [ x , . . . , x n ] . Then there exists a univariate polynomial P ( y ) ∈ R [ y ] and two vectors α ∈ N n , λ ∈ Z n \ { } such that p = x α P ( x λ ) if and only if for each pair ( i , j ) with ≤ i < j ≤ n, there is a polynomial P i j ( y ) ∈ R [ x , . . . , x i − , x i + , . . . , x j − , x j + , . . . , x n ][ y ] andfour integers β i j , β ji , µ i j , µ ji with µ i j , µ ji not both zero such that p = x β ij i x β ji j P i j ( x µ ij i x µ ji j ) .Proof. The necessity is clear. For the su ffi ciency, we proceed by induction on the number n ofvariables. There is nothing to show in the base case where n =
1. Assume that n > n − p as a polynomial in x , . . . , x n − over R [ x n ]. By the induction hypothesis, there isa polynomial P ∗ ( y ) ∈ R [ x n ][ y ] and two vectors ( α ∗ , . . . , α ∗ n − ) ∈ N n − , ( λ ∗ , . . . , λ ∗ n − ) ∈ Z n − withthe λ ∗ i not all zero such that p ( x n )( x , . . . , x n − ) = x α ∗ · · · x α ∗ n − n − P ∗ ( x λ ∗ · · · x λ ∗ n − n − ) . We may assume without loss of generality that λ ∗ ,
0. Regarding P ∗ ( y ) as an element of R [ y , x n ],we rewrite the preceding equation as p ( x , . . . , x n ) = x α ∗ · · · x α ∗ n − n − P ∗ ( x λ ∗ · · · x λ ∗ n − n − , x n ) . (5.1)By taking i = j = n in the assumption, we know that p = x β n x β n n P n ( x µ n x µ n n ) for P n ( y ) ∈ R [ x , . . . , x n − ][ y ] and β n , β n , µ n , µ n ∈ Z with µ n , µ n not both zero. Therefore, p ( q µ n x , x , . . . , x n − , q − µ n x n ) = cp ( x , . . . , x n ) with c = q β n µ n − β n µ n ∈ R . It then follows from (5.1) that P ∗ ( q µ n λ ∗ x λ ∗ · · · x λ ∗ n − n − , q − µ n x n ) = cq − µ n α ∗ P ∗ ( x λ ∗ · · · x λ ∗ n − n − , x n ), thatis, P ∗ ( q µ n λ ∗ y , q − µ n x n ) = cq − µ n α ∗ P ∗ ( y , x n ) . Applying Lemma 5.2 to P ∗ ( y , x n ) yields that there is a univariate polynomial P ( y ) ∈ R [ y ] and fourintegers α n , α ∗ n , λ n , λ ∗ n with λ n , λ ∗ n not both zero such that P ∗ ( y , x n ) = y α ∗ n x α n n P ( y λ ∗ n x λ n n ). Substituting y = x λ ∗ · · · x λ ∗ n − n − into this equation, together with (5.1), implies that p = x α P ( x λ ) with α = ( α ∗ + λ ∗ α ∗ n , . . . , α ∗ n − + λ ∗ n − α ∗ n , α n ) and λ = ( λ ∗ λ ∗ n , . . . , λ ∗ n − λ ∗ n , λ n ). The proof is concluded bynoticing that λ , .Inspired by the above proposition, we propose an algorithm which takes a multivariate poly-nomial as input and computes its q -integer linear decomposition in an iterative fashion so that at21ach iteration step, only two variables enter the game whereas the others are treated as coe ffi cientparameters. MultivariateQILD . Given a polynomial p ∈ R [ x , . . . , x n ], compute its q -integer linear decom-position.1. If p ∈ R then set c = p ; and return c .2. Set c = cont( p ) and f = prim( p ). If supp( f ) is a singleton then set α to be the only elementand update c = c f / x α ; and return c x α .3. If n = α to be the lowest degree of f with respect to x , m = λ m = P m ( y ) = f ( y ) / y α ; and return c x α Q mi = P i ( x λ i ).4. If n = BivariateQILD with input f ∈ R [ x , x ] to compute its q -integer linear decomposition f = x α x α P m Y i = P i ( x λ i x λ i );and then return c x α x α P Q mi = P i ( x λ i x λ i ).5. Set α = , P = m = g = cont x , x ( f ), and update f = prim x , x ( f ).6. If g , g ∈ R [ x , . . . , x n ], returning g = x ˜ α · · · x ˜ α n n ˜ P m Y i = ˜ P i ( x ˜ λ i · · · x ˜ λ in n ) , update α = α + (0 , , ˜ α , . . . , ˜ α n ), P = P ˜ P , and for i = , . . . , ˜ m iteratively update m = m + λ m = (0 , , ˜ λ i , . . . , ˜ λ in ), P m ( y ) = ˜ P i ( y ).7. If supp( f ) is a singleton then set α ∗ to be the only element and update α = α + α ∗ , c = c f / x α ∗ ; and return c x α P Q mi = P i ( x λ i ).8. Set Λ = { (cid:0) (1) , f ( y , x , . . . , x n ) (cid:1) } .For k = , . . . , n − Λ k + = {} .8.2 For (cid:0) ( µ , . . . , µ k ) , h ( y , x k + , . . . , x n ) (cid:1) in Λ k doCall the algorithm BivariateQILD with input h ∈ R [ x k + , . . . , x n ][ y , x k + ]to compute its q -integer linear decomposition h = y α ∗ x β ∗ k + P ∗ m ∗ Y i = P ∗ i ( y λ ∗ i x µ ∗ i k + , x k + , . . . , x n ) , (5.2)where P ∗ ∈ R [ y , x k + , . . . , x n ] and P ∗ i ( y , x k + , . . . , x n ) ∈ R [ y , x k + , . . . , x n ];then update α by adding the vector ( µ α ∗ , . . . , µ k α ∗ , β ∗ , , . . . , P by multiplying P ∗ ( x µ · · · x µ k k , x k + , . . . , x n ) and update Λ k + by joining theelements (cid:0) ( µ λ ∗ i , . . . , µ k λ ∗ i , µ ∗ i ) , P ∗ i ( y , x k + , . . . , x n ) (cid:1) for i = , . . . , m .22. Set g ∈ R [ x , . . . , x n ] to be the denominator of P . Update P to be its numerator, α i = α i − deg x i ( g ) for i = , . . . , n −
1, and for (cid:0) µ , h ( y ) (cid:1) in Λ n iteratively update m = m + λ m = µ and P m ( y ) = h ( y ).10. Return c x α P Q mi = P i ( x λ i ). Theorem 5.4.
Let p ∈ R [ x , . . . , x n ] . Then the algorithm MultivariateQILD correctly com-putes the q-integer linear decomposition of p.Proof. The correctness immediately follows from Proposition 5.3.
Theorem 5.5.
Let p ∈ Z [ q , q − ][ x , . . . , x n ] . Assume that the numerator and denominator of phave maximum degree d in each variable from { q , x , . . . , x n } separately and let || p || ∞ = β . Thenthe algorithm MultivariateQILD computes the q-integer linear decomposition of p over Z using O ∼ ( d n + + d n + log β ) word operations with classical arithmetic and O ∼ ( d n + + d n + log β ) with fast arithmetic.Proof. Let T ( n , d , log β ) denote the number of word operations used by the algorithm applied to p . The first three steps are exactly the same as the first approach introduced in the precedingsection. Thus we know that Step 1 takes no word operations, Step 2 uses O ∼ ( d n + + d n + log β )word operations with classical arithmetic and O ∼ ( d n + + d n + log β ) with fast arithmetic, andStep 3 gives that the initial cost T (1 , d , log β ) is in O ∼ ( d + d log β ) with classical arithmetic andO ∼ ( d + d log β ) with fast arithmetic. Step 4 deals with the bivariate case. By Theorem 3.11, thisstep yields that T (2 , d , log β ) is in O( d + d log β ) with classical arithmetic and O ∼ ( d + d log β )with fast arithmetic.In Step 5, by Lemma 3.10, the computation of the content and primitive part can be donewithin O ∼ ( d n + + d n + log β ) word operations with classical arithmetic and O ∼ ( d n + + d n + log β )with fast arithmetic. Notice that g ∈ Z [ q , x , . . . , x n ] has maximum degree at most d in eachvariable separately and max-norm of word length O( nd + log β ) by Lemma 3.9. Then Step 6takes O( T ( n − , d , nd + log β )) word operations. Step 7 takes linear time in the cardinalityof supp( f ), which is at most ( d + n . In Step 8, notice that for the k th iteration, the polyno-mial h ∈ Z [ q , x k + , . . . , x n ][ y , x k + ] has maximum degree at most d in each variable separatelyand max-norm of word length O( nd + log β ). Thus by Theorem 3.11 with v = n − k −
1, the k th iteration requires O ∼ ( d n − k + + d n − k + log β ) word operations with classical arithmetic andO ∼ ( d n − k + + d n − k + log β ) with fast arithmetic. Since 1 ≤ k ≤ n −
1, this step in total takesO ∼ ( d n + + d n + log β ) word operations with classical arithmetic and O ∼ ( d n + + d n + log β ) withfast arithmetic, dominating the costs of Steps 9 and 10.In summary, we obtain the recurrence relationO( T ( n , d , log β )) ⊂ O ∼ ( d n + + d n + log β ) + O( T ( n − , d , nd + log β )) , along with T (1 , d , log β ) ∈ O( d + d log β ) and T (2 , d , log β ) ∈ O( d + d log β ) with classicalarithmetic or O( T ( n , d , log β )) ⊂ O ∼ ( d n + + d n + log β ) + O( T ( n − , d , nd + log β )) , along with T (1 , d , log β ) ∈ O ∼ ( d + d log β ) and T (2 , d , log β ) ∈ O ∼ ( d + d log β ) with fastarithmetic. The announced cost follows. 23 O x •••••••••••••• (0 ,
15) (30 , z = x − x O x ••••••••••••••••••••• (0 , ,
14) (22 , z = x x − x O x ••••••••• (0 ,
4) (14 ,
0) (29 , Figure 3: Newton polygons constructed in the three stages in Example 5.6.
Example 5.6.
Consider the same polynomial p given by (4.1) as Example 4.8. In order tocompute its q-integer linear decomposition over Z [ q , q − ] , the algorithm MultivariateQILD (mainly Step 8) proceeds in the following three stages with their respective Newton polygonsplotted in Figure 3. Firstly, by viewing p as a polynomial in x , x over Z [ q , q − , x , x ] , applyingthe algorithm BivariateQILD to p givesp = x P (1) ( x − x , x , x ) (5.3) withP (1) ( y , x , x ) = qy x x + qy x + q y x x + qy x x + qy x x + q y x x − y x x − y x x − qy x x + q y x x − y x x + q y x x − y x x + q y x x − qy x x + qy x + qy x x + q y x − q y x x + qy x x − q y x x + qy x x − q y x x + q y x x + q yx x + q yx x + q x x . There is only one q-integer linear type, namely ( − , , of p over Z [ q , q − , x , x ] . Next, with inputP (1) ( y , x , x ) ∈ Z [ q , q − , x ][ y , x ] , calling the algorithm BivariateQILD again and substitutingy = x − x yields p = x · P · P (2) ( x x − x , x ) , (5.4) where P = qx x + x x + x x and P (2) ( y , x ) = q y x − q y x + qy x + qy + q y x − y x − y x + qyx + qx . The vector (2 , − , is then the only q-integerlinear type of p over Z [ q , q − , x ] . Finally, the last call to the algorithm BivariateQILD withinput P (2) ( y , x ) ∈ Z [ q , q − ][ y , x ] , along with the substitution y = x x − x , leads to the desireddecomposition (4.2) . The two q-integer linear types (2 , − , , and ( − , − , − , of p over Z [ q , q − ] have been correctly recovered.From (5.3) and (5.4) , one sees that p is q-integer linear over Z [ q , q − , x , x ] but it is notq-integer linear over Z [ q , q − , x ] . The latter in turn indicates the non-q-integer linearity of pover Z [ q , q − ] , even before starting the third stage. Once more, similar to the bivariate algorithm, the above algorithm can be modified so as todetermine the q -integer linearity of a given polynomial only. In other words, the algorithm canhalt already and return the negative answer whenever one of the following situations occurs.24 In Step 4 or in any iteration step of Step 8.2, any of the triggers of the bivariate algorithmlisted in Remark 3.8 is touched. • In Step 6, the polynomial g turns out to be not q -integer linear.
6. Complexity comparison
In this section, we discuss two more algorithms for computing the q -integer linear decompo-sition of polynomials, along with their complexity analyses in the case of bivariate polynomialsover Z [ q , q − ], so as to compare with our algorithms presented in Sections 4 and 5.The first algorithm is based on resultants and was developed by Le (2001), which itself servesas a q -analogue of the algorithm of Abramov and Le (2002) in the ordinary shift case. As alreadymentioned in the introduction, this algorithm is completely focused on bivariate polynomials. Sowe will further extend it to also tackle polynomials having more than two variables. The secondalgorithm is based on full irreducible factorization and work for polynomials in any numberof variables. This algorithm can be viewed as a q -analogue of the algorithm of Li and Zhang(2013) from the ordinary shift case. In order to give complexity comparison, we need to analyzethe costs of these two algorithms. As such we will briefly describe their main ideas. As we proceed with our bivariate algorithm, the algorithm of Le (2001) first finds candidatesfor q -integer linear types of a given bivariate polynomial and then obtains the corresponding uni-variate polynomials by going through these candidates. The di ff erence is that they use resultantsto determine candidates and perform bivariate GCD computations iteratively for detecting eachcandidate.In order to state its main idea, let p ∈ R [ x , y ] be a nonconstant polynomial which is primitivewith respect to its either variable. Then we know that p admits the q -integer linear decompositionof the form (3.1), in which all the λ i , µ i are nonzero. By Lemma 5.2, an integer pair ( λ, µ ) with λµ , q -integer linear types ( λ i , µ i ) if and only if there exists a nonconstant factor f ∈ R [ x , y ] of p with the property that f divides f ( q µ x , q − λ y ) in R [ x , y ]. Note that such an f must satisfy deg x ( f ) deg y ( f ) > f ( x , f (0 , y ) , p is assumed to be primitivewith respect to its either variable. By a careful study on the structure of the factor f , it is then nothard to see that f divides f ( q µ x , q − λ y ) in R [ x , y ] if and only if f divides f ( qx , q − λ/µ y ) in R [ x , y ].Observe that any integer pair ( λ, µ ) with λµ , r = − λ/µ .We have thus shown the following. Lemma 6.1.
With p given above, a nonzero rational number r gives rise to a q-integer lineartype of p if and only if gcd( p , p ( qx , q r y )) is nonconstant. This implies that all rationals r i = − λ i /µ i must be roots of the resultant Res y ( p , p ( qx , q r y )) ∈ R [ q r , x ] in terms of r , or equivalently, they are eliminated by the content in R [ q r ] of the resul-tant with respect to x . Note that such rational roots of a polynomial in R [ q r ] can be found bymatching powers of q appearing in the given polynomial in pairs along with a subsequent sub-stitution for zero testing. One can find more details in (Le, 2001, §
5) or (Le et al., 2001, § r i (and then for the ( λ i , µ i )). Aftergenerating candidates, the algorithm of Le (2001) continues to find the corresponding univariatepolynomials by calculating a factor f of p that stabilizes gcd( f , f ( qx , q r y )), or more e ffi ciently,25cd( f , f ( q µ x , q − λ y )) for each candidate r = − λ/µ . This operation actually induces bivariate poly-nomial arithmetic over R and thus may take considerably more time than Step 6.1 of our bivariatealgorithm BivariateQILD . In order to improve the performance, we instead proceed by usingStep 6 of our bivariate algorithm.Note that Lemma 6.1 cannot be literally carried over to a polynomial in more than two vari-ables. It is not clear how to directly generalize the algorithm of Le (2001) to the multivariate case.Nevertheless, as indicated by Proposition 5.3, this algorithm extends to the case of polynomialsin any number of variables in the same fashion as our bivariate algorithm.The lemma below provides bounds for the resultant of two trivariate integer polynomials,which can be verified by following the proof of (Bistritz and Lifshitz, 2010, Theorem 10) butarguing from the perspective of trivariate polynomials.
Lemma 6.2.
Let f , g ∈ Z [ q , x , y ] with deg q ( f ) = d q , deg x ( f ) = d x , deg y ( f ) = d y , deg q ( g ) = e q , deg x ( g ) = e x and deg y ( g ) = e y . Let R = Res y ( f , g ) . Then R ∈ Z [ q , x ] with deg q ( R ) ≤ d q e y + d y e q , deg x ( R ) ≤ d x e y + d y e x and || R || ∞ ≤ ( d y + e y )!(max { d q , e q } + d y + e y − (max { d x , e x } + d y + e y − || f || e y ∞ || g || d y ∞ . The following theorem gives a complexity analysis for the algorithm of Le (2001) whenapplied to a polynomial in Z [ q , q − ][ x , y ]. Theorem 6.3.
Let p be a polynomial in Z [ q , q − ][ x , y ] whose numerator and denominator havemaximum degree d in each variable from { q , x , y } and with || p || ∞ = β . Then the algorithm ofLe takes O ∼ ( d + d log β + d log β ) word operations with classical arithmetic and O ∼ ( d + d log β ) with fast arithmetic.Proof. With a slight abuse of notation, let p be the input polynomial with content with respectto its either variable being removed. Then p ∈ Z [ q , x , y ] and log || p || ∞ ∈ O( d + log β ). The algo-rithm proceeds to compute the resultant Res y ( p , p ( qx , q r y )) with r undetermined. Observe thatevery entry in the Sylvester matrix is a monomial in q r . Thus we have || Res y ( p , p ( qx , q r y )) || ∞ ≤|| Res y ( p , p ( qx , y )) || ∞ . By Lemma 6.2, Res y ( p , p ( qx , q r y )) has degree in q at most 3 d , degree in q r at most d , degree in x at most 2 d and max-norm at most B = (2 d )!(2 d + d − ( d + d − || p || d ∞ .Then log B ∈ O( d + d log β ). Viewing q r as a new indeterminate u independent of q , we can com-pute this resultant using a small prime modular algorithm, along with an evaluation-interpolationscheme: (1) choose ⌈ log (2 B + ⌉ primes, each of size O(log B ); (2) for every chosen prime h , dothe following: reduce all coe ffi cients of P ( x , y ) and P ( qx , uy ) modulo h , evaluate both modularimages at D = d points for ( q , u , x ), compute D univariate resultants of two polynomials in Z h [ y ] of degrees in y at most d , and recover the modular resultant by interpolation; (3) recon-struct the desired resultant using the Chinese remainder theorem. Ignoring the cost for choosingprimes in Step (1), we analyze the costs used by Steps (2)-(3). In Step (2), the cost per prime h forreducing all coe ffi cients modulo h is O( d log β log h ) word operations with classical arithmetic.The evaluation and interpolation steps are performed in O( d D ) ⊂ O( d ) arithmetic operations in Z h with classical arithmetic and O ∼ ( dD ) ⊂ O ∼ ( d ) with fast arithmetic. Each univariate resultantover Z h [ y ] can be computed in O( d ) arithmetic operations in Z h with classical arithmetic andO ∼ ( d ) with fast arithmetic, yielding O( d ) arithmetic operations in Z h with classical arithmeticand O ∼ ( d ) with fast arithmetic in total for this step. Notice that the cost for each arithmeticoperations in Z h is O(log h ) word operations with classical arithmetic and O ∼ (log h ) with fastarithmetic. Also notice that every chosen prime h is of word length log h ∈ O(log log B ). Thus26tep (2) in total takes O ∼ ( d log B ) word operations with classical arithmetic and O ∼ ( d log B )with fast arithmetic. In Step (3), the Chinese remainder theorem in total requires O ∼ ( d log B )word operations with classical arithmetic and O ∼ ( d log B ) with fast arithmetic. Therefore, com-puting the resultant Res y ( p , p ( qx , q r y )) takes O ∼ ( d + d log β + d log β ) word operations withclassical arithmetic and O ∼ ( d + d log β ) with fast arithmetic. This will dominate the costs forsubsequent steps including finding the rational roots and computing corresponding univariatepolynomials. The claimed cost follows. Similar to the ordinary shift case (Li and Zhang, 2013), the q -integer linear decomposition ofa multivariate polynomial can also be computed by full irreducible factorization. The key obser-vation is that, for any q -integer linear polynomial p ∈ R [ x , . . . , x n ] of only one type ( λ , . . . , λ n ),the di ff erence of any two vectors from supp( p ) can be written into the form k · ( λ , . . . , λ n ) forsome k ∈ Z . This allows one to readily determine the q -integer linearity of any irreduciblepolynomial. That is, given an irreducible polynomial p ∈ R [ x , . . . , x n ], take α ∈ N n to be theminimal vector of supp( p ) and investigate whether the di ff erence between α and any other vec-tor from supp( p ) is equal to a scalar multiple of the same integer vector. One thus immediatelyestablishes a factorization-based algorithm for the computation of the q -integer linear decom-position of a polynomial in R [ x , . . . , x n ]: first perform the full irreducible factorization of theinput polynomial over R ; then determine the q -integer linearity of each irreducible factor; finallyregroup all factors of the same q -integer linear type.A careful study of the above algorithm leads to the following complexity. Theorem 6.4.
Let p be a polynomial in Z [ q , q − ][ x , y ] whose numerator and denominator havemaximum degree d in each variable from { q , x , y } and with || p || ∞ = β . Then the factorization-based algorithm described above requires O ∼ ( d log β ) word operations with classical arith-metic and O ∼ ( d log β ) with fast arithmetic.Proof. Computing a complete factorization of p into irreducibles over Z [ q , q − ] dominates theother costs of the algorithm. This is essentially the complexity of factoring in Z [ q ][ x , y ], forpolynomials bounded by degree d in all variables ( q , x and y ). While we do not know of anexplicit analysis of this complexity (beyond being in polynomial-time, since (Kaltofen, 1985)),the algorithm of Gao (2003) can be applied and analyzed over the function field Q ( q ), and ap-pears to require O ∼ ( d log β ) word operations using classical arithmetic, and O ∼ ( d log β ) wordoperations using fast arithmetic.
7. Implementation and timings
We have implemented all algorithms in M aple R is thering of polynomials over Z [ q , q − ]. The code is available by email request. In order to get an ideaabout the e ffi ciency of our algorithms, we have compared their runtime, as well as the memoryrequirements, to the performance of our Maple implementations of the two algorithms discussedin the preceding section.The test suite was generated by p = P m Y i = num( P i ( x λ i )) , (7.1)where n , m ∈ N , 27 P ∈ Z [ q ][ x , . . . , x n ] is a random polynomial with deg x ,..., x n ( P ) = deg q ( P ) = d , • the λ i ∈ Z n are random integer vectors each of which has entries of maximum absolutevalue no more than 10 (note that they may not be distinct), • P i ( z ) = f i ( z ) f i ( z ) with f i j ( z ) ∈ Z [ q ][ z ] a random polynomial of degree j · d for some d ∈ N , and num( · · · ) denotes the numerator of the argument.Note that, in all tests, the algorithms take the expanded forms of examples given above as input.All timings are measured in seconds on a Linux computer with 128GB RAM and fifteen 1.2GHzDual core processors. The computations for the experiments did not use any parallelism.For a selection of random polynomials of the form (7.1) for di ff erent choices of n , m , d , d ,Table 1 collects the timings of the algorithm of Le (LQILD), the algorithm based on factorization(FQILD) and our two algorithms (MQILD , MQILD ). The dash in the table indicates that withthis choice of ( m , n , d , d ), the corresponding procedure reached the CPU time limit (which wasset to 12 hours) and yet did not return.( n , m , d , d ) LQILD FQILD MQILD MQILD (2 , , ,
1) 5408.48 0.04 0.01 0.01(2 , , ,
1) 8381.99 0.06 0.03 0.03(2 , , ,
1) – 0.19 0.04 0.04(2 , , ,
1) – 0.63 0.09 0.09(2 , , ,
1) – 1.47 0.13 0.10(2 , , ,
1) – 2.55 0.24 0.21(2 , , ,
1) – 6.64 0.42 0.39(2 , , ,
1) – 0.92 0.10 0.08(2 , , ,
1) – 3.29 0.31 0.26(2 , , ,
1) – 5.74 0.67 0.54(2 , , ,
1) – 18.83 2.01 1.54(2 , , ,
2) – 4.55 0.27 0.20(2 , , ,
2) – 114.82 4.98 4.53(2 , , ,
2) – 264.02 25.63 24.29(2 , , ,
2) – 36.14 1.38 1.21(2 , , ,
3) – 169.13 4.28 3.80(2 , , ,
4) – 649.03 12.15 12.86(2 , , ,
5) – 1554.31 31.54 33.50(2 , , ,
1) – 0.32 0.05 0.05(3 , , ,
1) – 1.99 0.14 0.12(4 , , ,
1) – 11.46 0.35 0.20(5 , , ,
1) – 183.17 0.99 0.63(6 , , ,
1) – 1141.32 2.58 0.98(7 , , ,
1) – 11759.89 6.07 1.74(8 , , ,
1) – 18153.45 10.60 5.29(9 , , ,
1) – – 65.53 38.12(10 , , ,
1) – – 176.25 89.87
Table 1:
Comparison of all four algorithms for a collection of polynomials p of the form (7.1). . Conclusion In this paper we have presented two new algorithms for computing the q -integer linear de-composition of a multivariate polynomial over any UFD of characteristic zero. When restrictedto the bivariate case, both algorithms reduce to the same bivariate algorithm. For the sake ofcomparison, we included an algorithm based on full irreducible factorization of polynomials.Compared with the known algorithm of Le (2001) and this factorization-based algorithm in thebivariate case, our algorithm is considerably faster. In practice, both our algorithms are also moree ffi cient than these two algorithms. In addition, we have extended and improved the original con-tribution of Le and provided complexity analysis for the improved version. We remark that bothour algorithms have much better performances than the other two algorithms in the case wherethe coe ffi cient domain contains algebraic numbers. References