Algorithms for Linearly Recurrent Sequences of Truncated Polynomials
AAlgorithms for Linearly Recurrent Sequencesof Truncated Polynomials
Seung Gyu Hyun
University of Waterloo
Waterloo, ON, Canada
Vincent Neiger
Univ. Limoges, CNRS, XLIM, UMR 7252
F-87000 Limoges, France
รric Schost
University of Waterloo
Waterloo, ON, Canada
ABSTRACT
Linear recurrent sequences are those whose elements are definedas linear combinations of preceding elements, and finding relationsof recurrences is a fundamental problem in computer algebra. Inthis paper, we focus on sequences whose elements are vectors overthe ring A = K [ ๐ฅ ]/โจ ๐ฅ ๐ โฉ of truncated polynomials. We presentthree methods for finding the ideal of canceling polynomials: aBerlekamp-Massey-like approach due to Kurakin, one which com-putes the kernel of some block-Hankel matrix over A via a minimalapproximant basis, and one based on bivariate Padรฉ approximation.We propose complexity improvements for the first two methods,respectively by avoiding the computation of redundant sequencegenerators and by exploiting the Hankel structure to compress theapproximation instance. We then observe these improvements inour empirical results through a C++ implementation. Finally wediscuss applications to the computation of minimal polynomialsand determinants of sparse matrices over A . KEYWORDS
Linear recurrences; Berlekamp-Massey-Sakata; Approximant basis;Sparse matrix.
Linear recurrences appear in many domains of computer scienceand mathematics, and computing recurrence relations efficiently is afundamental problem in computer algebra. More specifically, givena sequence of elements in K ๐ for some field K and positive integer ๐ , we seek a representation of its annihilator , which is a polynomialideal corresponding to all recurrence relations which are satisfied bythe sequence; the polynomials in the annihilator are said to cancel the sequence. In dimension ๐ = , the well-known Berlekamp-Massey algorithm [3, 27] computes the unique monic univariatepolynomial of minimal degree that cancels the sequence. Sakatafirst extended the Berlekamp-Massey algorithm to dimension [35] and then to the multi-dimensional case in general, i.e. with ๐ > [36]; see also the multi-dimensional extension by Norton andFitzpatrick [13]. More recent work includes variants of Sakataโsalgorithm such as one which handles relations that satisfy severalsequences simultaneously [37], approaches relating the problem tothe kernel of a multi-Hankel matrix and exploiting either fast linearalgebra [4] or a process similar to Gram-Schmidt orthogonalization[28], and an algorithm relying directly on multivariate polynomialarithmetic [5]. As for the representation of the annihilator, all thesealgorithms compute a Grรถbner basis or a border basis of the idealof recurrence relations.In this paper, we focus on computing recurrence relations forsequences whose elements are in A ๐ , where A = K [ ๐ฅ ]/โจ ๐ฅ ๐ โฉ . This problem can be solved using a specialization of Kurakinโsalgorithm [21, 22], as detailed in Section 3, where we also explicitlydescribe the output generating set of the annihilator as a lexico-graphic Grรถbner basis of some bivariate ideal. We derive a costbound of ๐ ห ( ๐ฟ๐ ( ๐ ๐ฟ๐ + ๐ ๐ ๐ )) operations in the base field K , where ๐ฟ is the order of the recurrence, and ๐ is an exponent for matrixmultiplication over K [8, 25]. Because the output Grรถbner basis isoften non-minimal, in Section 4 we modify Kurakinโs algorithm toavoid as much as possible the computation of redundant genera-tors, which lowers the cost to ๐ ห ( ๐ฟ๐ โ ( ๐ ๐ฟ๐ + ๐ ๐ ๐ )) , where ๐ โ is anumber arising in the algorithm as an estimation for the cardinality ๐ opt of minimal Grรถbner bases of the annihilator. In Section 7, weobserve that empirically ๐ โ is indeed often close or equal to ๐ opt .Despite the improvement, the above cost bound still has a de-pendence at least quadratic in the dimension ๐ . Our interest in thecase ๐ โซ is motivated among others by the following fact: givena zero-dimensional ideal I โ K [ ๐ฅ, ๐ฆ ] , one can recover a Grรถbnerbasis of it via I = Ann ( ๐ ) for some well-chosen ๐ โ A N only if K [ ๐ฅ, ๐ฆ ]/I has the Gorenstein property [16, 26]. When that is notthe case, one can recover a basis of I via the annihilator of several sequences simultaneously, which means precisely ๐ > . For large ๐ , we compute the annihilator via a minimal approximant basisof a block-Hankel matrix over A constructed from ๐ . Computingthis approximant basis via the algorithm PM-Basis of [14] leads toa complexity of ๐ ห ( ๐ฟ ๐ ๐๐ ) operations in K (Section 5.1). We thenpropose a novel improvement of this minimal approximant basiscomputation, based on a randomized compression of the input ma-trix which leverages its block-Hankel structure, reducing the costto ๐ ห ( ๐ฟ ๐๐ + ๐ฟ ๐ ๐ ) operations in K (Section 5.2).The four above algorithms have been implemented in C++ usingthe libraries NTL and PML, using Lazardโs structural theorem forgenerating examples of sequences; see Section 7 for more details.Our experiments on a prime field K highlight a good match betweencost bounds and practical running times, confirming also the benefitobtained from the improvements of both Kurakinโs algorithm andthe plain approximant basis approach.Furthermore, in Section 6 we propose an algorithm with costquasi-linear in the order ๐ฟ , whereas the above cost bounds are atleast quadratic. For ๐ โ ๐ ( ๐ฟ ) , we compute the annihilator via thebivariate Padรฉ approximation algorithm of [29]: this uses ๐ ห ( ๐ ๐ + ๐ฟ ) operations in K , at the price of restricting to ๐ โ ๐ ( ) .Finally, in Section 8 we mention applications to the computationof minimal polynomials and determinants of sparse matrices over A . To design Wiedemann-like algorithms [41] for such matrices ๐ด โ A ๐ ร ๐ , we need to compute annihilators from sequences ofthe form ( ๐ข ๐ ๐ด ๐ ๐ฃ ) ๐ โฅ โ A N for some vectors ๐ข and ๐ฃ ; several suchsequences may be needed, leading to the case ๐ > . a r X i v : . [ c s . S C ] F e b akataโs -dimensional algorithm shares similarities with thecase ๐ = of Kurakinโs algorithm, and has the same complexity ๐ ( ๐ฟ ๐ ) [35, Thm. 3]. Apart from this, to the best of our knowledgeprevious work has ๐ = and considers ๐ -dimensional sequencesover K for an arbitrary ๐ โฅ [4, 5, 28]. Complexity in this ๐ -variatecontext is often expressed using the degree ๐ท of the consideredzero-dimensional ideal; here, ๐ฟ โค ๐ท โค ๐ฟ๐ and a minimal Grรถbnerbasis or a border basis will have at most min ( ๐ฟ, ๐ ) + elements.The Scalar-FGLM algorithm has cost ๐ ห ( ๐ opt ๐ฟ ๐ ๐ ) [4, Prop. 16].Both the Artinian border basis and Polynomial-Scalar-FGLM algo-rithms [5, 28] cost ๐ ( ๐ท ๐ฟ๐ ) , which is ๐ ( ๐ฟ ๐ ) in the most favourablecase ๐ท = ๐ฟ , and ๐ ( ๐ฟ ๐ ) when ๐ท โ ฮ ( ๐ฟ๐ ) (which will be the casein our experiments, see Section 7). In all cases, a better complexitybound can be achieved by one of our algorithms outlined above.While this is not reflected in the cost estimates above, Kurakinโsalgorithm and our modified version are still affected by the shapeof the staircase of the computed Grรถbner basis, due to early termi-nation of the iterations and late additions; we leave a more refinedcomplexity analysis with respect to ๐ท as future work. In this section, we review key facts about linearly recurrent se-quences and algorithmic tools used throughout the paper. K [ ๐ฅ ]/โจ ๐ฅ ๐ โฉ We consider the set S = ( A ๐ ) N of (vector) sequences over thering A = K [ ๐ฅ ]/โจ ๐ฅ ๐ โฉ for some ๐ โ Z > , that is, sequences ๐ = ( ๐ , ๐ , . . . ) with each ๐ ๐ in A ๐ . Such a sequence is said to be lin-early recurrent if there exist ๐พ โ N and ๐ , . . . , ๐ ๐พ โ A with ๐ ๐พ invertible such that ๐ ๐ ๐ + ยท ยท ยท + ๐ ๐พ โ ๐ ๐ + ๐พ โ + ๐ ๐พ ๐ ๐ + ๐พ = for all ๐ โฅ (1)the order of ๐ is the smallest such ๐พ , denoted by ๐ฟ hereafter. Apolynomial ๐ + ยท ยท ยท + ๐ ๐พ ๐ฆ ๐พ in A [ ๐ฆ ] is said to cancel ๐ if ๐ , . . . , ๐ ๐พ satisfies Eq. (1) (without requiring that ๐ ๐พ be invertible). The set ofcanceling polynomials forms an ideal Ann ( ๐ ) in A [ ๐ฆ ] , called the annihilator of ๐ . Thus ๐ is linearly recurrent of order ๐ฟ if and onlyif there is a monic polynomial of degree ๐ฟ in Ann ( ๐ ) : such polyno-mials are called generating polynomials of ๐ . Unlike for sequencesover fields, here there may be canceling polynomials of degree lessthan ๐ฟ , which prevents uniqueness of generating polynomials; andthere are sequences which are not linearly recurrent but still admita nonzero canceling polynomial (i.e. Ann ( ๐ ) โ { } ). Example 2.1.
Consider A = K [ ๐ฅ ]/โจ ๐ฅ โฉ and the sequence ๐ = ( , + ๐ฅ, , + ๐ฅ, , + ๐ฅ, . . . ) in A N . Note that ๐ฅ ๐ = ( ๐ฅ, ๐ฅ, ๐ฅ, ๐ฅ, . . . ) .This sequence has order ๐ฟ = , a generating polynomial is ๐ฆ โ ,and a canceling polynomial of degree less than is ๐ฅ ( ๐ฆ โ ) . One canverify that Ann ( ๐ ) = โจ ๐ฆ โ , ๐ฅ ( ๐ฆ โ )โฉ ; in particular ๐ฆ + ๐ฅ ( ๐ฆ โ ) โ is also a generating polynomial. For any sequence ๐ in K N whichis not linearly recurrent, the sequence ๐ฅ ๐ in A N is not linearlyrecurrent but is cancelled by ๐ฅ , i.e. ๐ฅ โ Ann ( ๐ ) \ { } .Canceling polynomials can be characterized as denominatorsof the (vector) generating series of the sequence, defined as ๐บ = (cid:205) ๐ โฅ ๐ ๐ ๐ฆ โ ๐ โ in ( A [[ ๐ฆ โ ]]) ๐ . In what follows, the elements of A [ ๐ฆ ] ๐ are called polynomials, and for ๐ = ( ๐ , . . . , ๐ ๐ ) โ A [ ๐ฆ ] ๐ wedefine deg ( ๐ ) = max โค ๐ โค ๐ deg ( ๐ ๐ ) . Lemma 2.2. Let ๐ โ S , let ๐บ be its generating series, and let ๐ โ A [ ๐ฆ ] . Then, ๐ โ Ann ( ๐ ) if and only if the series ๐๐บ โ ( A [[ ๐ฆ โ ]]) ๐ is a polynomial, in which case deg ( ๐๐บ ) < deg ( ๐ ) . Proof. One has ๐๐บ = (cid:205) โค ๐ โค ๐พ (cid:205) ๐ โฅโ ๐ ๐ ๐ ๐ ๐ + ๐ ๐ฆ โ ๐ โ , where ๐พ = deg ( ๐ ) and ๐ = ๐ + ยท ยท ยท + ๐ ๐พ ๐ฆ ๐พ . Thus all terms of ๐๐บ have degreeless than ๐พ , and ๐๐บ is a polynomial if and only if its term in ๐ฆ โ ๐ โ vanishes for all ๐ โฅ , i.e. if and only if Eq. (1) holds. โก In this paper, we want to compute a generating set for
Ann ( ๐ ) ,but we typically only have access to a finite number of terms ofthe sequence for algorithms. Suppose we have access to the partialsequence ๐ ๐ = ( ๐ , . . . , ๐ ๐ โ ) in S ๐ = ( A ๐ ) ๐ , for some ๐ โ Z > .Similar to Eq. (1), a polynomial ๐ + ยท ยท ยท + ๐ ๐พ ๐ฆ ๐พ of degree ๐พ < ๐ cancels ๐ ๐ if ๐ ๐ ๐ + ยท ยท ยท + ๐ ๐พ ๐ ๐ + ๐พ = for all โค ๐ < ๐ โ ๐พ. (2)The next lemma shows that polynomials of degree ๐พ which cancel ๐ ๐ also cancel the whole sequence ๐ , provided the discrepancy between ๐ and ๐พ is sufficiently large (namely, ๐ โฅ ๐พ + ๐ฟ ).Lemma 2.3. Let ๐ โ S be linearly recurrent of order ๐ฟ . For any ๐ โ Z > and any ๐ โ A [ ๐ฆ ] with deg ( ๐ ) โค ๐ โ ๐ฟ , one has ๐ โ Ann ( ๐ ) if and only if ๐ cancels ๐ ๐ . Proof. Obviously, any polynomial ๐ โ Ann ( ๐ ) also cancels ๐ ๐ ,for any ๐ โ Z > greater than the degree of ๐ . Now let ๐ โ A [ ๐ฆ ] \{ } such that ๐พ = deg ( ๐ ) โค ๐ โ ๐ฟ and ๐ cancels ๐ ๐ . Since ๐ โ ๐พ โฅ ๐ฟ ,Eq. (2) yields (cid:205) โค ๐ โค ๐พ ๐ ๐ ๐ ๐ + ๐ = for โค ๐ < ๐ฟ . Furthermore, since ๐ is linearly recurrent of order ๐ฟ , there exists ๐ฆ ๐ฟ โ (cid:205) โค ๐ < ๐ฟ โ ๐ ๐ ๐ฆ ๐ โ A [ ๐ฆ ] which cancels ๐ , meaning that ๐ ๐ + ๐ = (cid:205) โค ๐ < ๐ฟ ๐ ๐ ๐ ๐ โ ๐ฟ + ๐ + ๐ forany ๐ โฅ ๐ฟ . Therefore we get โ๏ธ โค ๐ โค ๐พ ๐ ๐ ๐ ๐ + ๐ = โ๏ธ โค ๐ < ๐ฟ ๐ ๐ โ๏ธ โค ๐ โค ๐พ ๐ ๐ ๐ ๐ โ ๐ฟ + ๐ + ๐ . Using this identity, it follows by induction on ๐ โฅ ๐ฟ that the relation (cid:205) โค ๐ โค ๐พ ๐ ๐ ๐ ๐ + ๐ = also holds for all ๐ โฅ ๐ฟ . Hence ๐ โ Ann ( ๐ ) . โก Uni-dimensional sequences of vectors in A ๐ as above can be inter-preted as two-dimensional sequences of vectors in K ๐ , that is, se-quences ๐ = ( ๐ ๐,๐ ) ๐,๐ โฅ in ๐ = ( K ๐ ) N . This is based on the naturalinjective morphism ๐ : A [ ๐ฆ ] โ K [ ๐ผ, ๐ฝ ] with ( ๐ ( ๐ฅ ) , ๐ ( ๐ฆ )) = ( ๐ผ, ๐ฝ ) .Here we recall from [13, 35] that a polynomial ๐ = (cid:205) ๐,๐ ๐ ๐ ๐ ๐ผ ๐ ๐ฝ ๐ in K [ ๐ผ, ๐ฝ ] is said to cancel a sequence ๐ = ( ๐ ๐,๐ ) ๐,๐ โฅ โ ๐ if โ๏ธ ๐,๐ โฅ ๐ ๐ ๐ ๐ ๐ + ๐ ,๐ + ๐ = for all ๐ , ๐ โฅ . Then let ๐ = ( ๐ , ๐ , . . . ) โ S , and define ๐ = ( ๐ ๐,๐ ) ๐,๐ โฅ โ ๐ suchthat ๐ ๐,๐ โ K ๐ is the coefficient of degree ๐ โ โ ๐ of the truncatedpolynomial vector ๐ ๐ โ A ๐ if ๐ < ๐ , and ๐ ๐,๐ = otherwise. Then apolynomial ๐ โ A [ ๐ฆ ] cancels ๐ if and only if the polynomial ๐ ( ๐ ) cancels ๐ . Furthermore, ๐ is linearly recurrent if and only if theset of polynomials in K [ ๐ผ, ๐ฝ ] which cancel ๐ is a zero-dimensionalideal of K [ ๐ผ, ๐ฝ ] which contains ๐ผ ๐ .In what follows, we define ยฏ ๐ (I) = โจ{ ๐ ( ๐ ) | ๐ โ I} โช { ๐ผ ๐ }โฉ forany ideal I of A [ ๐ฆ ] , providing a correspondence between the ideals f A [ ๐ฆ ] and those of K [ ๐ผ, ๐ฝ ] containing ๐ผ ๐ . For insight into possibleโniceโ generating sets for Ann ( ๐ ) , we consider the lexicographicorder โผ lex with ๐ผ โผ lex ๐ฝ , and use the fact that Grรถbner bases ofthe ideals in K [ ๐ผ, ๐ฝ ] for this order are well understood [24]. Below,unless mentioned otherwise, we use โผ lex when some term order isneeded, e.g. leading terms and Grรถbner bases.Consider a zero-dimensional ideal I in K [ ๐ผ, ๐ฝ ] that contains apower of ๐ผ and let G be its reduced Grรถbner basis. Let ( ๐ฝ ๐ , ๐ผ ๐ ๐ฝ ๐ , . . . , ๐ผ ๐ ๐ก โ ๐ฝ ๐ ๐ก โ , ๐ผ ๐ ๐ก ) be the leading terms of the elements of G listed in decreasingorder, i.e. the ๐ ๐ โs are decreasing and the ๐ ๐ โs are increasing. Weset ๐ = ๐ ๐ก = , and for โค ๐ โค ๐ก we set ๐ฟ ๐ = ๐ ๐ โ ๐ ๐ โ , so that ๐ ๐ = ๐ฟ + ยท ยท ยท + ๐ฟ ๐ . Similarly, for โค ๐ < ๐ก we set ๐ ๐ = ๐ ๐ โ ๐ ๐ + .Then write G = { ๐ , . . . , ๐ ๐ก } , with ๐ ๐ having leading term ๐ผ ๐ ๐ ๐ฝ ๐ ๐ ;in particular ๐ ๐ก = ๐ผ ๐ ๐ก = ๐ผ ๐ฟ +ยทยทยท+ ๐ฟ ๐ก and ๐ is monic in ๐ฆ .Lazardโs Theorem states the following [24]: for โค ๐ โค ๐ก onecan write ๐ ๐ = ๐ผ ๐ ๐ ^ ๐ ๐ , with ^ ๐ ๐ monic of degree ๐ ๐ in ๐ฝ . In addition,for โค ๐ < ๐ก , ^ ๐ ๐ = ๐ ๐ / ๐ผ ๐ ๐ is in the ideal generated by โจ ^ ๐ ๐ + , ๐ผ ๐ฟ ๐ + ^ ๐ ๐ + , . . . , ๐ผ ๐ฟ ๐ + +ยทยทยท+ ๐ฟ ๐ก โฉ = (cid:28) ๐ ๐ + ๐ผ ๐ ๐ + , ๐ ๐ + ๐ผ ๐ ๐ + , . . . , ๐ ๐ก ๐ผ ๐ ๐ + (cid:29) ; in particular, ๐ผ ๐ฟ divides ๐ , . . . , ๐ ๐ก . Lazard also proved that a setof polynomials which satisfies these conditions is necessarily aminimal Grรถbner basis.With the above notation, a minimal Grรถbner basis of I has cardi-nality ๐ก + , with ๐ก โค min ( ๐ , ๐ ๐ก ) since = ๐ < ๐ < ยท ยท ยท < ๐ ๐ก and = ๐ ๐ก < ยท ยท ยท < ๐ < ๐ . Since for the reduced Grรถbner basis G eachpolynomial ๐ ๐ is represented by at most ๐ ๐ ๐ก coefficients in K , thetotal size of G in terms of field elements is at most ๐ ๐ ๐ก min ( ๐ , ๐ ๐ก ) .Finer bounds for the cardinality and size of G could be given usingthe vector space dimension dim K ( K [ ๐ผ, ๐ฝ ]/I) . For a univariate polynomial matrix ๐น โ K [ ๐ฅ ] ๐ ร ๐ and a positiveinteger ๐ , the set A ๐ ( ๐น ) = { ๐ โ K [ ๐ฅ ] ร ๐ | ๐๐น = ๐ฅ ๐ } is a free K [ ๐ฅ ] -module of rank ๐ whose elements are called approx-imants for ๐น at order ๐ [1, 40]. Bases of such submodules can berepresented as ๐ ร ๐ nonsingular matrices over K [ ๐ฅ ] and are usu-ally computed in specific forms, namely reduced forms [20, 42] andtheir canonical form called the Popov form [20, 33]. Extensions ofthese forms have been defined to accommodate degree weights ordegree constraints, and are called shifted reduced or Popov forms[1, 2, 40]. The algorithm PM-Basis [14] computes an approximantbasis in shifted reduced form using ๐ ห ( ๐ ๐ โ ( ๐ + ๐ ) ๐ ) operations in K ; using essentially two calls to this algorithm, one can recover theunique approximant basis in shifted Popov form within the samecost bound [19].More generally, in the bivariate case with ๐น โ K [ ๐ผ, ๐ฝ ] ๐ ร ๐ and ( ๐, ๐ ) โ Z > , the set A ๐,๐ ( ๐น ) = { ๐ โ K [ ๐ผ, ๐ฝ ] ร ๐ | ๐๐น = ๐ฅ ๐ , ๐ฆ ๐ } is a K [ ๐ผ, ๐ฝ ] -submodule of K [ ๐ผ, ๐ฝ ] ร ๐ whose elements are called approximants for ๐น at order ( ๐, ๐ ) [12, 32]. Such submodules areusually represented by a โผ -Grรถbner basis for some term order โผ on K [ ๐ผ, ๐ฝ ] ร ๐ ; for definitions of term orders and Grรถbner bases forsubmodules we refer to [9, 11]. For ๐ โค ๐ algorithms based on aniterative approach or on efficient linear algebra yield cost boundsin ๐ ห ( ๐ ( ๐๐๐ ) + ( ๐๐๐ ) ) and ๐ ห ( ๐ ( ๐๐๐ ) ๐ โ + ( ๐๐๐ ) ๐ ) operationsin K respectively [12, 31], whereas a recent divide and conquerapproach costs ๐ ห (( ๐ ๐ + ๐ ๐ ) ๐๐ ) field where ๐ = ๐ min ( ๐, ๐ ) [29, Prop. 5.5]; in these cases the output is a minimal Grรถbner basis. In [21], Kurakin gives an algorithm based on the Berlekamp-Masseyalgorithm that computes the annihilators of a partial sequence overa ring ๐ (and modules over ๐ ) that can be decomposed as a disjointunion ๐ = { } โช ๐ โช ยท ยท ยท โช ๐ ๐ โ where ๐ ๐ = { ๐ ๐ ๐ โ | ๐ โ โ ๐ invertible } for some ๐ ๐ โ ๐ . In this paper we consider ๐ = A = K [ ๐ฅ ]/โจ ๐ฅ ๐ โฉ ; in this case thecanonical choice is ๐ ๐ = ๐ฅ ๐ , with ๐ ๐ = { ๐ฅ ๐ ๐ โ | ๐ โ โ A with nonzero constant term } . Given a partial sequence ๐ ๐ of length ๐ over ( A ๐ ) ๐ , Kurakinโsalgorithm computes ๐ polynomials ๐ ๐ โ A [ ๐ฆ ] , ๐ = , . . . , ๐ โ , suchthat ๐ ๐ is a canceling polynomial of ๐ ๐ that has leading coefficient ๐ฅ ๐ and is minimal in degree among all canceling polynomials with lead-ing coefficient ๐ฅ ๐ . Furthermore, one has Ann ( ๐ ) = โจ ๐ , . . . , ๐ ๐ โ โฉ provided ๐ โฅ ๐ฟ [22, Thm. 1].We first define three operations on sequences. Given a partialsequence ๐ ๐ and ๐ โ A , ๐ ยท ๐ ๐ denotes multiplying ๐ to every elementin ๐ ๐ , while ๐ฆ ๐ ยท ๐ ๐ denotes a shift of ๐ elements โ that is, removingthe first ๐ elements. Given another partial sequence ^ ๐ ^ ๐ , the sum ๐ ๐ + ^ ๐ ^ ๐ returns the first min ( ๐, ^ ๐ ) elements of the two sequencesadded together element-wise.Kurakinโs algorithm iterates on ๐ = , . . . , ๐ โ , keeping trackof polynomials ๐ ๐,๐ and partial sequences ๐ ๐,๐,๐ , such that ๐ ๐,๐,๐ = ๐ ๐,๐ ยท ๐ ๐ = (cid:205) ๐ โ ๐ ๐ = ๐ ( ๐ ) ๐,๐ ยท ๐ฆ ๐ ยท ๐ ๐,๐,๐ , where ๐ ( ๐ ) ๐,๐ is the ๐ -th coefficient of ๐ ๐,๐ . We also have the invariant that the leading coefficient of ๐ ๐,๐ is ๐ฅ ๐ for all ๐ . For each index ๐ = , . . . , ๐ โ , the algorithm essentiallyattempts to either create a zero by using the partial sequences fromprevious iterations with equal number of leading zeros (similar toGaussian elimination), or shift the sequence if we cannot cancelthis element.At each iteration ๐ , let I [ ๐ ] be the A -submodule of A ๐ generatedby the elements ๐ ๐,๐,๐ โฒ [ ๐ ] for all ๐ = , . . . , ๐ โ and ๐ โฒ < ๐ suchthat ๐ ๐,๐,๐ โฒ has ๐ leading zeros, and let P [ ๐, ๐ ] and S[ ๐, ๐ ] be thecorresponding polynomial and partial sequence to the ๐ -th elementin the basis of I [ ๐ ] . At iteration ๐ , if ๐ ๐,๐,๐ has ๐ leading zeros and ๐ ๐,๐,๐ [ ๐ ] โ I [ ๐ ] , then we can find coefficients such that ๐ ๐,๐,๐ [ ๐ ] โ (cid:205) ๐ ๐ ๐ I [ ๐, ๐ ] = and ๐ ๐,๐,๐ โ (cid:205) ๐ ๐ ๐ S[ ๐, ๐ ] results in a sequencewith ๐ + zeros since both sequences had ๐ leading zeros and wecanceled ๐ ๐,๐,๐ [ ๐ ] . The algorithm terminates when all ๐ ๐,๐,๐ = (seeAlgorithm 1).We track the subiterations by the index ๐ก for analysis; this doesnot play a role in the algorithm. Kurakin shows that the total num-ber of subiterations across all ๐ is ๐ ( ๐ ) per polynomial, bringing thetotal to ๐ ( ๐๐ ) ([21, Thm. 2]). However, the analysis of the runtime lgorithm 1 Kurakinโs Algorithm ( ๐ ๐ ) Input: partial sequence ๐ ๐ Output: minimal canceling polynomials of ๐ ๐ for ๐ = , . . . , ๐ โ do set ๐ ๐, = ๐ฅ ๐ and ๐ ๐,๐, = ๐ฅ ๐ ๐ ๐ set ๐ to be index of first non-zero element of ๐ ๐,๐, if ๐ ๐,๐, [ ๐ ] = then continue to next ๐ if ๐ ๐,๐, [ ๐ ] โ then add ๐ ๐,๐, [ ๐ ] , ๐ ๐, , ๐ ๐,๐, to I [ ๐ ] , P [ ๐ ] , S[ ๐ ] resp. for ๐ = , . . . , ๐ โ do for ๐ = , . . . ๐ โ do set ๐ก = set ๐ ( ๐ก ) ๐,๐ = ๐ฆ๐ ๐,๐ โ set ๐ ( ๐ก ) ๐,๐,๐ = ๐ฆ ยท ๐ ๐,๐,๐ โ (shift by 1) if ๐ ( ๐ก ) ๐,๐,๐ = then continue to next ๐ set ๐ to be the first non-zero index of ๐ ( ๐ก ) ๐,๐,๐ if ๐ ( ๐ก ) ๐,๐,๐ [ ๐ ] โ I [ ๐ ] then continue to next ๐ else solve for ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ] โ (cid:205) ๐ ๐ ๐ I [ ๐, ๐ ] = set ๐ ( ๐ก + ) ๐,๐,๐ = ๐ ( ๐ก ) ๐,๐,๐ โ (cid:205) ๐ ๐ ๐ S[ ๐, ๐ ] set ๐ ( ๐ก + ) ๐,๐ = ๐ ( ๐ก ) ๐,๐ โ (cid:205) ๐ ๐ ๐ P [ ๐, ๐ ] go to line 12 with ๐ก = ๐ก + for ๐ = , . . . , ๐ โ do set ๐ ๐,๐,๐ = ๐ ( ๐ก ) ๐,๐,๐ and ๐ ๐,๐ = ๐ ( ๐ก ) ๐,๐ set ๐ to be the index of first non-zero element of ๐ ๐,๐,๐ if ๐ ๐,๐,๐ [ ๐ ] โ I [ ๐ ] then add ๐ ๐,๐,๐ [ ๐ ] , ๐ ๐,๐ , ๐ ๐,๐,๐ to I [ ๐ ] , P [ ๐ ] , S[ ๐ ] resp. reduce the basis of I [ ๐ ] if needed for ๐ = , . . . , ๐ โ do return ๐ ๐,๐ that makes ๐ ๐,๐,๐ = for the first timein [21] treats all ring operations (including computing solutionto line 17 of Algorithm 1) as constants, which is unrealistic over A ๐ . Thus, we will give a cost analysis in terms of number of fieldoperations over K .We note that, since A ๐ is a free K [ ๐ฅ ] -module of rank ๐ (witha basis given by the canonical vectors of length ๐ ) and K [ ๐ฅ ] is aprincipal ideal domain, any of its K [ ๐ฅ ] -submodule is free of rank atmost ๐ . As a consequence, the number of generators of I [ ๐ ] is atmost ๐ . This will allow us to bound the cost for solving submodulemembership as well as the equation ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ] โ (cid:205) ๐ ๐ ๐ I [ ๐, ๐ ] = .We can check membership ๐ ๐,๐,๐ [ ๐ ] โ I [ ๐ ] and solve ๐ ๐,๐ ,๐ [ ๐ ] โ (cid:205) ๐ ๐ I [ ๐, ๐ ] = by finding the right approximant basis of ๐น = (cid:2) I [ ๐, ] ยท ยท ยท I [ ๐, ๐ โ ] ๐ ๐,๐ ,๐ [ ๐ ] (cid:3) in Popov form. Since ๐น has ๐ rows and at most ๐ + columns, wecan compute this in cost ๐ ห ( ๐ ๐ ๐ ) [19]. For lines 18 and 19, ๐ [ ๐, ๐ ] and ๐ [ ๐, ๐ ] have length and degree at most ๐ resp., making the costof these two lines ๐ ห ( ๐ ( ๐๐๐ )) = ๐ ห ( ๐ ๐๐ ) . Finally, using the factthat the total number of subiterations is bounded by ๐ ( ๐๐ ) , wearrive at the total cost of ๐ ห ( ๐๐ ( ๐ ๐๐ + ๐ ๐ ๐ )) operations over K . We conclude by showing that the output of Algorithm 1 is indeeda basis of Ann ( ๐ ) and that it forms a Grรถbner basis wrt lexicograph-ical order when viewed as bivariate polynomials.Theorem 3.1. For each ๐ โ { , . . . , ๐ โ } , let ๐ ๐ be a cancelingpolynomial of ๐ with leading coefficient ๐ฅ ๐ that is minimal in degreeamong all polynomials with leading coefficient ๐ฅ ๐ . Then one has Ann ( ๐ ) = โจ ๐ , . . . , ๐ ๐ โ โฉ . Furthermore, { ๐ ( ๐ ) , ยท ยท ยท , ๐ ( ๐ ๐ โ ) , ๐ผ ๐ } forms a Grรถbner basis of ยฏ ๐ ( Ann ( ๐ )) with respect to the lexicographicterm order with ๐ผ โผ lex ๐ฝ . Proof. Suppose that there exists some ๐ โ A [ ๐ฆ ] with leadingcoefficient ๐ฅ ๐ก that is in Ann ( ๐ ) but ๐ โ โจ ๐ , . . . , ๐ ๐ โ โฉ . Note thatfor any polynomial in A [ ๐ฆ ] , we can always make the leading coef-ficient to be some ๐ฅ ๐ก by pulling out the minimal power of ๐ฅ fromthe leading coefficient and multiplying by its inverse. Now, sincewe assumed minimality of degrees for ๐ ๐ โs, deg ( ๐ ) > deg ( ๐ ๐ก ) and ๐ โฒ = ๐ โ ๐ฆ deg ๐ โ deg ๐ ๐ก ๐ ๐ก โ Ann ( ๐ ) has degree less than ๐ . By nor-malizing the leading coefficient of ๐ โฒ to be some ๐ฅ ๐ก โฒ , we can repeatthe same process and keep decreasing the degree. This process mustterminate when we encounter some ๐ โฒ with leading coefficient ๐ฅ ๐ก โฒ such that deg ๐ โฒ < deg ๐ ๐ก โฒ , or ๐ โฒ = . Both cases lead to contradic-tions; thus, such ๐ cannot exist and Ann ( ๐ ) = โจ ๐ , . . . , ๐ ๐ โ โฉ .Next, let ๐บ = { ๐ , . . . , ๐ ๐ } , ๐ ๐ โ K [ ๐ผ, ๐ฝ ] with leading coeffi-cient ๐ฅ ๐ ๐ , be the minimal reduced (lexicographic) Grรถbner basisof ยฏ ๐ ( Ann ( ๐ )) . We can turn ๐บ into another non-minimal Grรถbnerbasis by adding the polynomials ๐ ๐ ๐ ๐ , for ๐ = , . . . , ๐ ๐ + โ ; wedefine the resulting basis as ๐บ โฒ = { ๐ โฒ , ยท ยท ยท , ๐ โฒ ๐ } , with ๐ โฒ ๐ = ๐ผ ๐ andeach ๐ โฒ ๐ has leading term ๐ผ ๐ ๐ฝ ๐ ๐ . Furthermore, define ๐ข ๐ as the degreeof ๐ ๐ such that ๐ ( ๐ ๐ ) has leading term ๐ผ ๐ ๐ฝ ๐ข ๐ .For ๐ = , . . . , ๐ , we have that ๐ข ๐ โฅ ๐ ๐ , otherwise ๐บ โฒ would notreduce ๐ ( ๐ ๐ ) to zero, which ๐บ โฒ must since ๐ ( ๐ ๐ ) โ ยฏ ๐ ( Ann ( ๐ )) . Wealso have that ๐ข ๐ โค ๐ ๐ due to the assumed minimality of degree for ๐ ๐ โs. Thus, the leading terms of { ๐ ( ๐ ) , . . . , ๐ ( ๐ ๐ โ ) , ๐ผ ๐ } generatethe leading terms of ยฏ ๐ ( Ann ( ๐ )) and forms a Grรถbner basis. โก Kurakinโs algorithm requires that we keep track of all ๐ possiblegenerators, regardless of the actual number of generators needed.For example, consider ๐ = ( , , , , , . . . ) โ A N . One can verifythat Ann ( ๐ ) = โจ ๐ฆ โ ๐ฆ โ โฉ , but Kurakinโs algorithm will return { ๐ฅ ๐ ( ๐ฆ โ ๐ฆ โ )} for all ๐ = , . . . , ๐ โ . In this section, we willoutline a modified version of Kurakinโs algorithm that attempts toavoid as many extraneous computations as possible.In the previous example, we can see that the polynomials asso-ciated with ๐ฅ ๐ , ๐ โฅ , were not useful. The next definition aims toqualify precisely the usefulness of the monomial ๐ฅ ๐ . Definition 4.1.
Let ๐ ๐,๐ and ๐ ๐,๐,๐ be the polynomial and sequenceat the end of step ๐ associated with monomial ๐ฅ ๐ . A monomial ๐ฅ ๐ is useful wrt to ๐ฅ ๐ , ๐ < ๐ , at step ๐ if at least one of two conditionsis true at the end of ๐ :U1. ๐ ๐ ,๐ โ ๐ฅ ๐ โ ๐ ๐ ๐ ,๐ U2. let ๐ ๐ and ๐ ๐ be the index of the first non-zero element of ๐ ๐,๐ ,๐ and ๐ ๐,๐ ,๐ resp., then ๐ ๐ โ ๐ ๐ uppose a monomial ๐ฅ ๐ is not useful wrt ๐ฅ ๐ at step ๐ , then bynegation condition U1, we have ๐ ๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ๐ ,๐ . Due to nega-tion of U2, ๐ ๐,๐ ,๐ is the zero sequence if and only if ๐ ๐,๐ ,๐ is thezero sequence; so either we return ๐ ๐ ,๐ = ๐ฅ ๐ก โ ๐ก ๐ ๐ ,๐ or we donot terminate at this step for both monomials. Finally, since ๐ ๐ = ๐ ๐ and ๐ ๐,๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ๐,๐ ,๐ , we always have that ๐ ๐,๐ ,๐ [ ๐ ๐ ] = ๐ฅ ๐ โ ๐ ๐ ๐,๐ ,๐ [ ๐ ๐ ] โ (cid:0) โจ ๐ ๐,๐ ,๐ [ ๐ ๐ ]โฉ โช I [ ๐ ๐ ] (cid:1) , meaning we can safelyignore ๐ ๐,๐ ,๐ [ ๐ ๐ ] when updating I [ ๐ ๐ ] at the end of step ๐ . Thus,the negation of usefulness conditions U1 and U2 implies that anycomputation associated with ๐ฅ ๐ is not needed at step ๐ .However, as defined, U1 and U2 do not impose any conditionsabout the subiterations (indexed by ๐ก ). The next lemma gives adifferent characterization of the usefulness conditions in terms ๐ก .Lemma 4.2. If ๐ฅ ๐ is useful wrt to ๐ฅ ๐ at some step ๐ , then at somesubiteration ๐ก of step ๐ , one of the following conditions is true at thestart of ๐ก : u1. ๐ ( ๐ก ) ๐ ,๐ โ ๐ฅ ๐ โ ๐ ๐ ( ๐ก ) ๐ ,๐ u2. if ๐ ( ๐ก ) ๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ( ๐ก ) ๐ ,๐ , then ๐ ( ๐ก ) ๐ โ ๐ ( ๐ก ) ๐ u3. if ๐ ( ๐ก ) ๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ( ๐ก ) ๐ ,๐ and ๐ ( ๐ก ) ๐ = ๐ ( ๐ก ) ๐ , then ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] and ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] Proof. Suppose the conditions u1-u2-u3 are all false for everysubiteration ๐ก at ๐ . The negation of u1 forces ๐ ( ๐ก ) ๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ( ๐ก ) ๐ ,๐ at the start of ๐ก , which sets the hypothesis of u2 true, implying ๐ ( ๐ก ) ๐ = ๐ ( ๐ก ) ๐ . Finally, since the hypothesis of u3 holds, we musthave ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] or ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] . The twoare mutually exclusive since ๐ ( ๐ก ) ๐,๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ( ๐ก ) ๐,๐ ,๐ , if ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โI [ ๐ ( ๐ก ) ๐ ] , then ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] . When ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] ,we can update ๐ ( ๐ก + ) ๐ ,๐ = ๐ ( ๐ก ) ๐ ,๐ โ โ๏ธ ๐ ๐ I [ ๐ ( ๐ก ) ๐ , ๐ ] ๐ ( ๐ก + ) ๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ( ๐ก ) ๐ ,๐ โ ๐ฅ ๐ โ ๐ โ๏ธ ๐ ๐ P [ ๐ ( ๐ก ) ๐ , ๐ ] = ๐ฅ ๐ โ ๐ ๐ ( ๐ก + ) ๐ ,๐ , which was already implied by the assumption that condition u1was false for all ๐ก .On the other hand, when ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] , we also have ๐ ( ๐ก ) ๐,๐ ,๐ [ ๐ ( ๐ก ) ๐ ] โ I [ ๐ ( ๐ก ) ๐ ] , so the subiterations terminate and we musthave ๐ ๐ ,๐ = ๐ฅ ๐ โ ๐ ๐ ๐ ,๐ with ๐ ๐ = ๐ ๐ . This implies U1 and U2 alsodo not hold for step ๐ . โก While the converse is not true, we say a monomial ๐ฅ ๐ is poten-tially useful wrt ๐ฅ ๐ when at some step ๐ and subiteration ๐ก , at leastone of the conditions u1-u3 holds. Rather than iterating through ๐ = , . . . , ๐ โ , we keep a list of potentially useful monomials U and iterate through ๐ โ U , with U = [ ] initially. At each subit-eration, we check to see if there exists ๐ โฒ > ๐, ๐ โฒ โ U such that ๐ฅ ๐ โฒ satisfies one of u2 or u3, and add the smallest such ๐ โฒ to U . Notethat we need not check u1 since if u1 holds, then either u2 or u3must have been true at some previous subiteration, thus ๐ โฒ is alreadyincluded in U . Condition u2 can be checked in ๐ ( ๐ ) by checkingthe valuations of all entries in ๐ ๐,๐,๐ [ ๐ ] at lines 5 and 13. Conditionu3 can be checked in ๐ ( log ๐ ) membership computations via abinary search to find the minimal ๐ โฒ such that ๐ฅ ๐ โฒ โ ๐ ๐ ๐,๐,๐ [ ๐ ] โ I [ ๐ ] when ๐ ๐,๐,๐ [ ๐ ] โ I [ ๐ ] on line 14. Thus, the complexity for the subit-erations do not change in terms of ๐ ห (ยท) . Defining ๐ โ = |U| โค ๐ ,this brings the total cost to ๐ ห ( ๐๐ โ ( ๐ ๐๐ + ๐ ๐ ๐ )) . While we do notknow how far ๐ โ is from the number of polynomials in the minimallexicographic Grรถbner basis of ยฏ ๐ ( Ann ( ๐ )) , ๐ opt , we have observedexperimentally that ๐ โ is often equal or close to ๐ opt (see Section7). Extending the classical theory of linearly recurrent sequences overthe field K , another approach is to consider the left kernel of theblock-Hankel matrix ๐ป ๐ ,๐ = ๏ฃฎ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฐ ๐ ๐ ยท ยท ยท ๐ ๐ โ ๐ ๐ . . . ๐ ๐ ... . . . . . . ...๐ ๐ ๐ ๐ + ยท ยท ยท ๐ ๐ โ ๏ฃน๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃป โ A ( ๐ + )ร( ๐๐ ) . Indeed, if ๐ is large enough, vectors in this kernel represent poly-nomials which cancel ๐ , and which even generate all of Ann ( ๐ ) .Lemma 5.1. Let ๐ โ S be linearly recurrent of order ๐ฟ , and define K ๐ ,๐ = { ๐ = ๐ + ยท ยท ยท + ๐ ๐ ๐ฆ ๐ โ A [ ๐ฆ ] | [ ๐ ยท ยท ยท ๐ ๐ ] ๐ป ๐ ,๐ = } for ๐ โ N . Assume ๐ โฅ ๐ฟ . Then K ๐ ,๐ = Ann ( ๐ ) โฉ A [ ๐ฆ ] โค ๐ , and inparticular K ๐ ,๐ is a generating set of Ann ( ๐ ) . Proof. Let ๐ = ๐ +ยท ยท ยท+ ๐ ๐ ๐ฆ ๐ โ A [ ๐ฆ ] and ๐พ = deg ( ๐ ) โค ๐ . Then ๐ โ K ๐ ,๐ if and only if [ ๐ ยท ยท ยท ๐ ๐ ] ๐ป ๐ ,๐ = , and by definition ofcanceling partial sequences this exactly means that ๐ cancels ๐ ๐ + ๐พ .Now, deg ( ๐ ) = ๐พ โค ๐ + ๐พ โ ๐ฟ holds under the assumption ๐ โฅ ๐ฟ ,hence ๐ cancels ๐ ๐ + ๐พ if and only if ๐ โ Ann ( ๐ ) by Lemma 2.3. Itfollows that K ๐ ,๐ generates Ann ( ๐ ) , since there exists a generatingset of Ann ( ๐ ) whose polynomials all have degree at most ๐ฟ . โก Computing the left kernel of ๐ป ๐ ,๐ can be done via univariate ap-proximation. Indeed, calling ๐น โ K [ ๐ฅ ] ( ๐ + )ร( ๐๐ ) the natural liftingof ๐ป ๐ ,๐ , an approximant basis of ๐น at order ๐ gives a generatingset of that left kernel; As recalled in Section 2.4, using PM-Basis, abasis of A ๐ ( ๐น ) in shifted reduced or Popov form can be computedin ๐ ห ( ๐ ๐ โ ( ๐ + ๐๐ ) ๐ ) = ๐ ห ( ๐ ๐ ๐๐ ) operations in K . Now we show that, when ๐ is large, one can speed-up the above ap-proach by a randomized โcompressionโ of the matrix ๐ป ๐ ,๐ . Precisely,taking a random constant matrix ๐ถ โ K ( ๐๐ )ร( ๐ + ) and performingthe right-multiplication ๐น๐ถ , one obtains a square ( ๐ + ) ร ( ๐ + ) matrix such that A ๐ ( ๐น ) = A ๐ ( ๐น๐ถ ) holds with good probability.The cost of the approximant basis computation is thus reduced to ๐ ห ( ๐ ๐ ๐ ) operations in K , and the right-multiplication can be doneefficiently by leveraging the block-Hankel structure of ๐น .Theorem 5.2. Algorithm 2 takes as input an integer ๐ โ Z > ,vectors ๐น , . . . , ๐น ๐ + ๐ โ โ K [ ๐ฅ ] ร ๐ of degree less than ๐ , and a shift ๐ค โ Z ๐ > , and uses ๐ ห ( ๐๐๐๐ + ๐ ๐ ๐ ) operations in K to compute a ๐ค -Popov matrix ๐ โ K [ ๐ฅ ] ๐ ร ๐ of degree at most ๐ . It chooses at most ๐๐๐ elements independently and uniformly at random from a subset f K of cardinality ๐ , and ๐ is the ๐ค -Popov basis of A ๐ ( ๐น ) withprobability at least โ ๐๐ , where ๐น is the block-Hankel matrix ๐น = ๏ฃฎ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฏ๏ฃฐ ๐น ๐น ยท ยท ยท ๐น ๐ โ ๐น ๐น . . . ๐น ๐ ... . . . . . . ...๐น ๐ โ ๐น ๐ ยท ยท ยท ๐น ๐ + ๐ โ ๏ฃน๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃบ๏ฃป โ K [ ๐ฅ ] ๐ ร( ๐๐ ) . (3)When applied to the computation of Ann ( ๐ ) with ๐ = ๐ + , thecost becomes ๐ ห ( ๐ ๐๐ + ๐ ๐ ๐ ) . Below we focus on the case of interest ๐ โค ๐๐ , since when ๐๐ โ ๐ ( ๐ ) this ๐ค -Popov approximant basis iscomputed deterministically by PM-Basis at a cost of ๐ ห ( ๐ ๐ ๐ ) oper-ations in K . Our approach is based on the following two lemmas.Lemma 5.3. Let ๐น โ K [ ๐ฅ ] ๐ ร ๐ and ๐ โ Z > . Let ๐ถ โ K [ ๐ฅ ] ๐ ร ๐ and ๐พ โ K [ ๐ฅ ] ๐ ร( ๐ โ ๐ ) , for some ๐ โ { , . . . , ๐ } , such that ๐น๐พ = and [ ๐ถ ( ) ๐พ ( )] โ K ๐ ร ๐ is invertible. Then, ๐ โฅ ๐ where ๐ is the rankof ๐น , and A ๐ ( ๐น ) = A ๐ ( ๐น๐ถ ) . Proof. Let ๐ = [ ๐ถ ๐พ ] โ K [ ๐ฅ ] ๐ ร ๐ . The assumption that ๐ ( ) is invertible ensures that ๐ is nonsingular (since det ( ๐ )( ) = det ( ๐ ( )) โ ), and therefore ๐พ has full rank ๐ โ ๐ . The assumptionthat the columns of ๐พ are in the right kernel of ๐น , which has rank ๐ โ ๐ , implies that ๐ โ ๐ โค ๐ โ ๐ and therefore ๐ โฅ ๐ .The inclusion A ๐ ( ๐น ) โ A ๐ ( ๐น๐ถ ) is obvious. For the other in-clusion, let ๐ โ A ๐ ( ๐น๐ถ ) , i.e. there exists ๐ โ K [ ๐ฅ ] ร ๐ such that ๐๐น๐ถ = ๐ฅ ๐ ๐ . It follows that ๐๐น ๐ = ๐ฅ ๐ [ ๐ ] , and thus ๐๐น = ๐ฅ ๐ [ ๐ ] ๐ โ = ๐ฅ ๐ [ ๐ ] Ajd ( ๐ ) det ( ๐ ) where Adj ( ๐ ) โ K [ ๐ฅ ] ๐ ร ๐ is the adjugate of ๐ . Our assumption det ( ๐ )( ) โ means that ๐ฅ ๐ and det ( ๐ ) are coprime, hence det ( ๐ ) divides [ ๐ ] Ajd ( ๐ ) , and ๐๐น = ๐ฅ ๐ follows. โก Lemma 5.4.
Let ๐น โ K [ ๐ฅ ] ๐ ร ๐ with rank ๐ and ๐ โค ๐ , and let ๐ โ{ ๐, . . . , ๐ } . Let R be a finite subset of K of cardinality ๐ โ Z > , and let ๐ถ โ K ๐ ร ๐ with entries chosen independently and uniformly at randomfrom R . Then, the probability that there exists ๐พ โ K [ ๐ฅ ] ๐ ร( ๐ โ ๐ ) suchthat [ ๐ถ ๐พ ( )] is invertible and ๐น๐พ = is at least โ ๐๐ ; furthermoreif K is finite and K = K , this probability is at least (cid:206) ๐๐ = ( โ ๐ โ ๐ ) . Proof. Consider a right kernel basis ๐ต โ K [ ๐ฅ ] ๐ ร( ๐ โ ๐ ) for ๐น .Then ๐ต has unimodular row bases [43, Lem. 3.1], implying thatthere exists ๐ โ K [ ๐ฅ ] ( ๐ โ ๐ )ร ๐ such that ๐ ๐ต = ๐ผ ๐ โ ๐ . In particular ๐ ( ) ๐ต ( ) = ๐ผ ๐ โ ๐ and therefore ๐ต ( ) has full rank ๐ โ ๐ . Define ๐พ โ K [ ๐ฅ ] ๐ ร( ๐ โ ๐ ) as the matrix formed by the first ๐ โ ๐ columns of ๐ต (recall ๐ โ ๐ โค ๐ โ ๐ by assumption). Then ๐น๐พ = . Furthermore ๐พ ( ) has rank ๐ โ ๐ , hence the DeMillo-Lipton-Schwartz-Zippellemma implies that [ ๐ถ ๐พ ( )] โ K ๐ ร ๐ is singular with probabilityat most ๐ / ๐ [10, 38, 44]. If K is finite and K = K then [ ๐ถ ๐พ ( )] isinvertible with probability exactly (cid:206) ๐๐ = ( โ ๐ โ ๐ ) . โก These lemmas lead to Algorithm 2 and Theorem 5.2; indeedcomputing ๐น๐ถ has quasi-linear cost ๐ ห ( ๐๐๐๐ ) thanks to the block-Hankel structure of ๐น , and then the call PM-Basis ( ๐, ๐น๐ถ, ๐ค ) costs ๐ ห ( ๐ ๐ ๐ ) operations as recalled in Section 2.4.Note that โ ๐ / ๐ โฅ / as soon as ๐ โฅ ๐ (which implies ๐ โฅ ๐ );furthermore (cid:206) ๐๐ = ( โ ๐ โ ๐ ) โฅ / already for ๐ = . The random-ization is of the Monte Carlo type, since the algorithm may return Algorithm 2
Hankel-PM-Basis ( ๐, ๐น, ๐ค ) Input: integers ๐, ๐, ๐, ๐ โ Z > , vectors ๐น , . . . , ๐น ๐ + ๐ โ โ K [ ๐ฅ ] ร ๐ of degree less than ๐ , a shift ๐ค โ Z ๐ > Output: a ๐ค -Popov matrix ๐ โ K [ ๐ฅ ] ๐ ร ๐ of degree at most ๐ ๐น โ K [ ๐ฅ ] ๐ ร( ๐๐ ) โ form the block-Hankel matrix as in Eq. (3) if ๐ โฅ ๐๐ then return PM-Basis ( ๐, ๐น, ๐ค ) Choose ๐ โ { ๐, . . . , ๐ } where ๐ is the rank of ๐น (by default,choose ๐ = ๐ if no information is known on ๐ ) Fill a matrix ๐ถ โ K ( ๐๐ )ร ๐ with entries chosen uniformly andindependently at random from a subset of K of cardinality ๐ Compute ๐น๐ถ โ K [ ๐ฅ ] ๐ ร ๐ (exploiting the Hankel structure of ๐น ) return PM-Basis ( ๐, ๐น๐ถ, ๐ค ) ๐ which is not a basis of A ๐ ( ๐น ) . Still, since the expected ๐ค -Popovbasis ๐ of A ๐ ( ๐น ) is unique, one can easily increase the probabilityof success by repeating the randomized computation and followinga majority rule. Another approach is to rely on the non-interactive,Monte Carlo certification protocol of [15], which has lower costthan Algorithm 2 but requires a larger field K ; this first asks tocompute the coefficient of degree ๐ of ๐๐น , which here can be donevia bivariate polynomial multiplication in time ๐ ห ( ๐๐๐๐ ) thanks tothe structure of ๐น . For a given output ๐ , this certification can berepeated for better confidence in ๐ (in which case the coefficient ofdegree ๐ of ๐๐น needs only be computed once). Now, we propose another approach which directly uses the inter-pretation of cancelling polynomials as denominators of the gen-erating series of the sequence (see Lemma 2.2). The next lemmadescribes more precisely the link between the annihilator and thesedenominators when we have access to a partial sequence, that is,denominators of the generating series truncated at some order.One can also view this lemma as a description of the kernel of theunivariate Hankel matrix ๐ป ๐ ,๐ via bivariate Padรฉ approximation.Lemma 6.1. Let ๐ โ S be linearly recurrent of order ๐ฟ , and for ๐ โ N define ๐บ = (cid:205) ๐ < ๐ ๐ ๐ ๐ฆ ๐ โ โ ๐ โ A [ ๐ฆ ] ๐ and P ๐ ,๐ = { ๐ โ A [ ๐ฆ ] โค ๐ | ๐๐บ = ๐ mod ๐ฆ ๐ for some ๐ โ A [ ๐ฆ ] ๐ < ๐ } . Assume ๐ โฅ ๐ฟ . Then P ๐ ,๐ = Ann ( ๐ ) โฉ A [ ๐ฆ ] โค ๐ , and in particular P ๐ ,๐ is a generating set of Ann ( ๐ ) ; furthermore for any ๐ โ P ๐ ,๐ thecorresponding ๐ โ A [ ๐ฆ ] ๐ < ๐ satisfies deg ( ๐ ) < deg ( ๐ ) . Proof. Let ๐ = ๐ + ยท ยท ยท + ๐ ๐พ ๐ฆ ๐พ โ A [ ๐ฆ ] โค ๐ where ๐พ = deg ( ๐ ) .Then ๐ โ P ๐ ,๐ if and only if the coefficient of ๐๐บ of degree ๐ โ โ ๐ are zero for โค ๐ < ๐ . Since ๐พ โค ๐ โค ๐ โ โ ๐ , this coefficient is Coeff ( ๐๐บ, ๐ โ โ ๐ ) = ๐พ โ๏ธ ๐ = ๐ ๐ ๐ ๐ โ โ( ๐ โ โ ๐ โ ๐ ) = ๐พ โ๏ธ ๐ = ๐ ๐ ๐ ๐ + ๐ = . Thus we have proved P ๐ ,๐ = K ๐ ,๐ , and Lemma 5.1 shows theclaims in this lemma except the last one. Let ๐ โ P ๐ ,๐ and de-fine ๐ as the polynomial in A [ ๐ฆ ] ๐ < ๐ such that ๐๐บ = ๐ mod ๐ฆ ๐ .Since ๐ โ Ann ( ๐ ) , Lemma 2.2 shows that ๐๐บ is a polynomial. Onthe other hand the definitions of ๐บ and ๐บ yield ๐๐บ = ๐ฆ ๐ ๐๐บ โ ๐ (cid:205) ๐ โฅ ๐ ๐ ๐ ๐ฆ ๐ โ โ ๐ . Hence โ ๐ (cid:205) ๐ โฅ ๐ ๐ ๐ ๐ฆ ๐ โ โ ๐ is a polynomial,and since it has degree less than ๐พ , and thus in particular less than ๐ , it is equal to ๐ . โก rom ๐บ , define ๐น โ K [ ๐ผ, ๐ฝ ] ร ๐ of bi-degree less than ( ๐, ๐ ) via the morphism ๐ from Section 2.3. Equip K [ ๐ผ, ๐ฝ ] with the lex-icographic order โผ lex , and let โผ be the corresponding term overposition order on K [ ๐ผ, ๐ฝ ] ๐ + . Then a minimal โผ -Grรถbner basis ofthe submodule of simultaneous Padรฉ approximants {( ๐, ๐ ) โ K [ ๐ผ, ๐ฝ ] ร K [ ๐ผ, ๐ฝ ] ร ๐ | ๐๐น = ๐ mod ๐ฅ ๐ , ๐ฆ ๐ } is computed in ๐ ห (( ๐ ๐ min ( ๐, ๐ ) ๐ + ๐ min ( ๐, ๐ ) ) ๐๐ ) operations,using the algorithm of [29] (see also Section 2.4) with input matrixof size ( ๐ + ) ร ๐ formed by stacking the identity ๐ผ ๐ below ๐น .Lemma 6.1 shows that from this โผ -Grรถbner basis one can find aminimal โผ lex -Grรถbner basis of ยฏ ๐ ( Ann ( ๐ )) by selecting ๐ for each ( ๐, ๐ ) in the basis such that deg ๐ฝ ( ๐ ) < deg ๐ฝ ( ๐ ) .While the PM-Basis approach had cost quasi-linear in ๐ and ๐ ,the method here is most efficient in an opposite parameter range:for ๐ โ ๐ ( ) and ๐ โค ๐ the above cost bound becomes ๐ ห ( ๐ ๐ + ๐ ) . In this section, we compare timings for our implementations of thealgorithms in Sections 3 to 5. The algorithms were implemented inC++ using the libraries NTL [39] and PML [18] which provide high-performance support for univariate polynomials and polynomialmatrices. We leave the implementation of the bivariate algorithmof Section 6 as future work.To control the cardinality and shape of the lexicographic Grรถb-ner basis, we use Lazardโs structural theorem recalled in Section 2.3.The shape of the monomial staircase is randomized with maximal ๐ฝ -degree ๐ฟ and ๐ผ ๐ included in the basis. After generating a randomGrรถbner basis ๐บ of target degree and size, we use it to generate ๐ sequences (with ๐ = ๐ฟ terms), using random initial conditions.Finally, we provide the sequence as input and compute the annihi-lator of the sequence, which may not necessarily recover ๐บ itself(see Section 8.1 for more details). Runtimes are showed in Table 1. ๐ฟ ๐ ๐ ๐ opt ๐ท / ๐๐ฟ K LK ๐ โ PM-B HPM256 64 1 1 1 62.8 0.93 1 1.06 NA256 64 1 39 0.43 26.8 1.31 42 2.14 NA256 64 1 49 0.62 38.0 1.65 53 2.10 NA512 128 1 16 0.92 >100 12 17 20.5 NA512 12 3 4 0.4 6.93 1.40 4 2.47 1.86256 17 2 2 0.5 14.1 0.91 2 0.33 0.29256 16 8 1 1 54.1 3.16 1 0.56 0.25256 16 32 1 1 >100 39.8 1 2.79 0.35128 16 64 1 1 >100 >100 1 1.02 0.1332 128 1 12 0.91 7.85 0.078 12 0.029 NA32 256 1 14 0.94 27.3 0.12 14 0.08 NA128 256 1 27 0.92 >100 1.28 27 1.60 NA256 512 1 29 0.96 >100 8.65 29 27.8 NA
Table 1:
Runtimes of algorithms Kurakin, Lazy Kurakin, direct PM-Basis, and Hankel-PM-Basis. All timings are taken on AMD Ryzen5 3600X 6-Core CPU with 16 GB RAM. The base field is K = F . As we claim in Section 4, we can see that ๐ โ is often closeor equal to ๐ opt . More interestingly, Lazy Kurakin outperformsKurakin more than ๐ / ๐ โ would suggest. For example, for ๐ฟ = , ๐ = , ๐ opt = , ๐ / ๐ โ โ but runtime of Kurakin is 20times slower than Lazy Kurakin. This is because the complexitybound of ๐ ห ( ๐๐ โ ( ๐ ๐๐ + ๐ ๐ ๐ )) for lazy Kurakin assumes that ๐ โ polynomials are tracked from the beginning of the algorithms. How-ever, due to its lazy nature, polynomials are often added later inthe algorithm and the bound of ๐๐ โ subiterations may significantlyoverestimate the true number of subiterations.When ๐ฟ, ๐, ๐ are fixed, Kurakinโs algorithm performs worse for ๐ opt = than ๐ opt > , although this is a favourable case for LazyKurakin. In this case, Kurakinโs algorithm computes ๐ ๐ = ๐ฅ ๐ ๐ sothere cannot be any early termination. Additionally, the size ofthe stair case is maximal ( ๐ท = ๐๐ ), so this is also the worst casefor algorithms that whose complexity depends directly on ๐ท . LazyKurakinโs algorithm somewhat remedies this by using the extrastructure of A and adding monomials in a lazy fashion. When itis known that Ann ( ๐ ) = โจ ๐ โฉ , it is possible to design an algorithmthat is quasilinear in ๐ via structured system solving.For scalar sequences over A (that is, the case ๐ = ), Lazy Ku-rakinโs algorithm seems to be the best choice when the order ๐ islarge compared to ๐ , whereas PM-Basis seems to be the best choicein the converse. When ๐ = ๐ฟ = ๐ , Lazy Kurakin outperformsPM-Basis, given that ๐ โ is small. This is predicted by the theoreti-cal complexities, as the former has complexity ๐ ห ( ๐ ๐ โ ) , while thelatter has complexity ๐ ห ( ๐ ๐ + ) .For ๐ > , PM-Basis and Hankel-PM-Basis clearly outperformKurakin and Lazy Kurakin. This is as predicted since the complexityof the former depends linearly on ๐ , while the latter has a factor ๐ ๐ .The theoretical improvement of Hankel-PM-Basis over PM-Basisis observed empirically, especially for the two cases of ๐ = , . In this section, we outline two applications to sparse matrices over A = K [ ๐ฅ ]/โจ ๐ฅ ๐ โฉ . Firstly, given a sparse matrix ๐ด โ A ๐ ร ๐ , we wantto compute its minimal polynomials, which are polynomials of min-imal degree that cancel the matrix sequence ๐ ๐ด = ( ๐ด , ๐ด , ๐ด , ยท ยท ยท ) .Secondly, given ๐ด as above, we want to compute its determinant.In what follows, we will assume ๐ด has sparsity ๐ ( ๐ ) (i.e. it has ๐ ( ๐ ) nonzero entries) and that the representation of ๐ด allows us tocompute matrix-vector products at cost ๐ ห ( ๐๐ ) . The algorithms arebased on the well-known Wiedemann algorithm [41], which com-putes generating polynomials for sequences over fields, modifiedto account for the added complexity of working over A . Given a matrix ๐ด , the well-known Cayley-Hamilton theorem statesthat ๐ด cancels its own characteristic polynomial. This implies thatthe sequence of successive powers of ๐ด is linearly recurrent, and apolynomial of minimal degree that cancels this sequence is said tobe a minimal polynomial of ๐ด . A different view one can take is thatsuch canceling polynomials must cancel the ๐ linearly recurrentsequences (( ๐ด ๐ ) ๐ ,๐ ) ๐ โฅ simultaneously for โค ๐ , ๐ โค ๐ . Then,as usual, we want compute a Grรถbner basis of the ideal of thesecanceling polynomials, denoted by Ann ( ๐ด ) .Over A , trying to deduce Ann ( ๐ด ) from Ann (( ๐ข ๐ ๐ด ๐ ๐ฃ ) ๐ โฅ ) , forrandom vectors ๐ข, ๐ฃ โ A ๐ ร , presents a problem when Ann ( ๐ด ) doesnot have the Gorenstein property [16, 26]. When
Ann ( ๐ด ) does havethe Gorenstein property, it has been showed that Ann ( ๐ด ) can berecovered, with high probability, by using a bidimensional sequencewith random initial conditions, given that the characteristic of K is arge [4]. When it does not, then Ann ( ๐ด ) is still recoverable with asimilar approach, but using several sequences [30]. Over variouscommutative rings, the problem of computing minimal polynomialsof a matrix have been studied in [7, 17, 34]. However, the algorithmsgiven in these works do not exploit sparsity.Given matrix ๐ด as above, we start by choosing random ๐ข , ๐ฃ โ A ๐ and generating ๐ ๐ด, = ( ๐ข ๐ ๐ด ๐ ๐ฃ ) โค ๐ < ๐ . Note that by the resultsin Section 2.3, ๐ ๐ด, can be viewed as a truncated bidimensionalsequence. Next, we apply one of the algorithms in the previoussections to compute Ann ( ๐ ๐ด, ) . If Ann ( ๐ ๐ด, ) = Ann ( ๐ด ) , whichcan be checked probabilistically by checking if Ann ( ๐ ๐ด, ) alsocancels some validation sequence (( ๐ข โฒ ) ๐ ๐ด ๐ ๐ฃ ) โค ๐ < ๐ , we terminatethe process. Otherwise, we double the number of sequences bydoubling the number of random ๐ข ๐ โs and generating ๐ ๐ด, , . . . , ๐ ๐ด, ๐ .The cost of the process is ๐ ห ( ๐๐ ๐ + L( ๐, ๐, ๐ )) , where ๐ is thenumber of sequences used and L( ๐, ๐, ๐ ) is the cost of finding theannihilators of a partial sequence of length ๐ in ( K [ ๐ฅ ]/โจ ๐ฅ ๐ โฉ) ๐ . Notethat this process must terminate. The crudest bound is when ๐ > ๐ since then we could simply compute Ann ( ๐ด ) directly. Anotherslighly more refined bound for number of generic linear formsneeded is ๐ โค ๐ท , where ๐ท is the size of the staircase of Ann ( ๐ด ) [30,Prop. 1]. Given a matrix, we can deduce its determinant from its minimalpolynomial only when the characteristic polynomial is equal tothe minimal polynomial. Wiedemann [41] calls such matrices non-derogatory and shows that preconditioning any matrix ๐ต โ K ๐ ร ๐ with a random diagonal matrix ๐ท results in a nonderogatory matrixwith high probability. We will show that the same preconditioningcan be applied to matrices over A .Here, a particular role will be played by sequences ๐ โ ( A ๐ ) N such that Ann ( ๐ ) = โจ ๐ โฉ , for some monic ๐ โ A [ ๐ฆ ] . Indeed, thenext theorem shows that it is sufficient for the constant part of ๐ด to be nonderogatory in K for ๐ด to be nonderogatory in A and forthe sequence of its powers to satisfy this property.Theorem 8.1. Let ๐ด โ K ๐ ร ๐ be the constant part of ๐ด (i.e. for ๐ฅ = ). If ๐ด is nonderogatory, then Ann ( ๐ด ) = โจ ๐ โฉ for some monic ๐ โ A [ ๐ฆ ] of degree ๐ . Proof. Let ๐ โ A [ ๐ฆ ] be the minimal monic polynomial of the se-quence ๐ ๐ด = ( ๐ด , ๐ด , ๐ด , . . . ) . Then deg ( ๐ ) โค ๐ , since consideringthe characteristic polynomial of ๐ด shows that the order of ๐ ๐ด is atmost ๐ . Since ๐ด is nonderogatory, any canceling polynomial musthave degree โฅ ๐ ; thus, deg ( ๐ ) = ๐ . Furthermore, we can show that ๐ ๐ = ๐ฅ ๐ ๐ is the unique minimal polynomial with leading coefficient ๐ฅ ๐ . If there exists another polynomial ๐ of degree ๐ and leadingcoefficient ๐ฅ ๐ such that ๐ โ ๐ฅ ๐ ๐ , then ๐ โ ๐ฅ ๐ ๐ is a canceling poly-nomial of degree less than ๐ , contradicting the previous statement.Thus, ๐, ๐ , . . . , ๐ ๐ โ are minimal in degree and, by Theorem 3.1, Ann ( ๐ด ) = โจ ๐, ๐ , . . . , ๐ ๐ โ โฉ = โจ ๐ โฉ . โก The above theorem allows us to use the same preconditioner asin [41]: a random constant diagonal matrix ๐ท . The preconditioningensures that the ideal of canceling polynomial is generated by asingle monic polynomial; thus, ยฏ ๐ ( Ann ( ๐ด๐ท )) is Gorenstein andrequires only a single linear form to be recovered. Furthermore, when it is known that the ideal is generated by a single polynomial,we can recover this polynomial in ๐ ห ( ๐๐ ) by taking advantage ofthe fact that the constant part of the leading ๐ ร ๐ submatrix of ๐ป ๐ , ๐ is an invertible Hankel matrix [6]. Once we have ๐ , we cancompute det ( ๐ด ) = ๐ ( )( (cid:206) ๐ ๐ท ๐,๐ ) โ .Under our sparsity assumption, the cost of this method is ๐ ห ( ๐ ๐ ) for computing ( ๐ข ๐ ๐ด ๐ ๐ฃ ) ๐ โค ๐ , ๐ ห ( ๐๐ ) for computing ๐ , and ๐ ห ( ๐ + ๐ ) for recovering the determinant from ๐ , leading to the total cost of ๐ ห ( ๐ ๐ ) operations in K . This is to be compared with computingthe determinant of ๐ด โat full precisionโ, i.e. by seeing ๐ด as a matrixover K [ ๐ฅ ] , and then truncating the result modulo ๐ฅ ๐ : this costs ๐ ห ( ๐ ๐ ๐ ) operations in K [23]. REFERENCES [1] B. Beckermann and G. Labahn. 1994. A Uniform Approach for the Fast Computa-tion of Matrix-Type Padรฉ Approximants.
SIAM J. Matrix Anal. Appl.
15, 3 (1994),804โ823. https://doi.org/10.1137/S0895479892230031[2] B. Beckermann, G. Labahn, and G. Villard. 1999. Shifted Normal Forms of Polyno-mial Matrices. In
ISSACโ99 . ACM, 189โ196. https://doi.org/10.1145/309831.309929[3] E. Berlekamp. 1968. Nonbinary BCH decoding (Abstr.).
IEEE Trans. Inf. Theory
14, 2 (1968), 242โ242. https://doi.org/10.1109/TIT.1968.1054109[4] J. Berthomieu, B. Boyer, and J.-C. Faugรจre. 2017. Linear algebra for computingGrรถbner bases of linear recursive multidimensional sequences.
J. Symb. Comput.
83 (2017), 36โ67. https://doi.org/10.1016/j.jsc.2016.11.005[5] J. Berthomieu and J.-C. Faugรจre. 2018. A Polynomial-Division-Based Algorithmfor Computing Linear Recurrence Relations. In
ISSACโ18 . 79โ86. https://doi.org/10.1145/3208976.3209017[6] A. Bostan, C.-P. Jeannerod, and ร Schost. 2008. Solving structured linear systemswith large displacement rank.
Theor. Comput. Sci.
Communications in Algebra
33, 12(2005), 4491โ4504. https://doi.org/10.1080/00927870500274820[8] D. Coppersmith and S. Winograd. 1990. Matrix multiplication via arithmeticprogressions.
J. Symb. Comput.
9, 3 (1990), 251โ280. https://doi.org/10.1016/S0747-7171(08)80013-2[9] D. A. Cox, J. Little, and D. OโShea. 2005.
Using Algebraic Geometry (second edition) .Springer-Verlag New-York, New York, NY. https://doi.org/10.1007/b138611[10] R. A. DeMillo and R. J. Lipton. 1978. A Probabilistic Remark on Algebraic ProgramTesting.
Inform. Process. Lett.
7, 4 (1978), 193โ195.[11] D. Eisenbud. 1995.
Commutative Algebra: with a View Toward Algebraic Geometry .Springer. https://doi.org/10.1007/978-1-4612-5350-1[12] P. Fitzpatrick. 1997. Solving a Multivariable Congruence by Change of Term Order.
J. Symb. Comput.
24, 5 (1997), 575โ589. https://doi.org/10.1006/jsco.1997.0153[13] P. Fitzpatrick and G. H. Norton. 1990. Finding a basis for the characteristic idealof an n-dimensional linear recurring sequence.
IEEE Trans. Inf. Theory
36, 6(1990), 1480โ1487. https://doi.org/10.1109/18.59953[14] P. Giorgi, C.-P. Jeannerod, and G. Villard. 2003. On the complexity of polynomialmatrix computations. In
ISSACโ03 . ACM, 135โ142. https://doi.org/10.1145/860854.860889[15] Pascal Giorgi and Vincent Neiger. 2018. Certification of Minimal ApproximantBases. In
ISSACโ18 . ACM, 167โ174. https://doi.org/10.1145/3208976.3208991[16] W. Grรถbner. 1935. รber irreduzible Ideale in kommutativen Ringen.
Math. Ann.
Linear Algebra Appl.
527 (2017), 12โ31. https://doi.org/10.1016/j.laa.2017.03.028[18] S. G. Hyun, V. Neiger, and ร. Schost. 2019. Implementations of Efficient Univari-ate Polynomial Matrix Algorithms and Application to Bivariate Resultants. In
ISSACโ19 . ACM, 235โ242. https://doi.org/10.1145/3326229.3326272[19] C.-P. Jeannerod, V. Neiger, and G. Villard. 2020. Fast computation of approximantbases in canonical form.
J. Symb. Comput.
98 (2020), 192โ224. https://doi.org/10.1016/j.jsc.2019.07.011[20] T. Kailath. 1980.
Linear Systems . Prentice-Hall.[21] V. L. Kurakin. 1998. The BerlekampโMassey algorithm over finite rings, modules,and bimodules.
Discrete Mathematics and Applications
8, 5 (1998), 441โ474.[22] V. L. Kurakin. 2000. Construction of the Annihilator of a Linear RecurringSequence over Finite Module with the help of the Berlekamp-Massey Algorithm.In
FPSAC 2000 . Springer, 476โ483. https://doi.org/10.1007/978-3-662-04166-6_45[23] G. Labahn, V. Neiger, and W. Zhou. 2017. Fast, deterministic computation of theHermite normal form and determinant of a polynomial matrix. 42 (2017), 44โ71.https://doi.org/10.1016/j.jco.2017.03.003
24] D. Lazard. 1985. Ideal Bases and Primary Decomposition: Case of Two Variables.
J. Symb. Comput.
1, 3 (1985), 261โ270.[25] F. Le Gall. 2014. Powers of Tensors and Fast Matrix Multiplication. In
ISSACโ14 (Kobe, Japan). ACM, 296โ303. https://doi.org/10.1145/2608628.2608664[26] F. S. Macaulay. 1934. Modern algebra and polynomial ideals. In
Math. Proc. Camb.Philos. Soc , Vol. 30. Cambridge University Press, 27โ46.[27] J. Massey. 1969. Shift-register synthesis and BCH decoding.
IEEE Trans. Inf.Theory
15 (1969), 122โ127.[28] B. Mourrain. 2017. Fast Algorithm for Border Bases of Artinian GorensteinAlgebras. In
ISSACโ17 (Kaiserslautern, Germany). ACM, 333โ340. https://doi.org/10.1145/3087604.3087632[29] S. Naldi and V. Neiger. 2020. A Divide-and-Conquer Algorithm for ComputingGrรถbner Bases of Syzygies in Finite Dimension. In
ISSACโ20 . ACM, 380โ387.https://doi.org/10.1145/3373207.3404059[30] V. Neiger, H. Rahkooy, and ร. Schost. 2017. Algorithms for zero-dimensionalideals using linear recurrent sequences. In
CASC 2017 . Springer, 313โ328.[31] V. Neiger and ร. Schost. 2020. Computing syzygies in finite dimension using fastlinear algebra.
J. Complexity
60 (2020), 101502. https://doi.org/10.1016/j.jco.2020.101502[32] H. OโKeeffe and P. Fitzpatrick. 2002. Grรถbner basis solutions of constrainedinterpolation problems.
Linear Algebra Appl.
351 (2002), 533โ551. https://doi.org/10.1016/S0024-3795(01)00509-2[33] V. M. Popov. 1972. Invariant Description of Linear, Time-Invariant ControllableSystems.
SIAM Journal on Control
10, 2 (1972), 252โ264. [34] R. Rissner. 2016. Null ideals of matrices over residue class rings of principal idealdomains.
Linear Algebra Appl.
494 (2016), 44โ69. https://doi.org/10.1016/j.laa.2016.01.004[35] S. Sakata. 1988. Finding a minimal set of linear recurring relations capable ofgenerating a given finite two-dimensional array.
J. Symb. Comput.
5, 3 (1988),321โ337. https://doi.org/10.1016/S0747-7171(88)80033-6[36] S. Sakata. 1990. Extension of the Berlekamp-Massey algorithm to N dimensions.
Information and Computation
84, 2 (1990), 207โ239.[37] S. Sakata. 2009. The BMS Algorithm. In
Grรถbner Bases, Coding, and Cryptography .Springer, 143โ163. https://doi.org/10.1007/978-3-540-93806-4_9[38] J. T. Schwartz. 1980. Fast Probabilistic Algorithms for Verification of PolynomialIdentities.
J. ACM
27, 4 (1980), 701โ717. https://doi.org/10.1145/322217.322225[39] V. Shoup. 2020. NTL: A Library for doing Number Theory, version 11.4.3. .[40] M. Van Barel and A. Bultheel. 1992. A general module theoretic framework forvector M-Padรฉ and matrix rational interpolation.
Numer. Algorithms
IEEETrans. Inf. Theory
32, 1 (1986), 54โ62. https://doi.org/10.1109/TIT.1986.1057137[42] W. A. Wolovich. 1974.
Linear Multivariable Systems . Applied MathematicalSciences, Vol. 11. Springer-Verlag New-York.[43] W. Zhou and G. Labahn. 2013. Computing Column Bases of Polynomial Matrices.In
ISSACโ13 . ACM, 379โ386. https://doi.org/10.1145/2465506.2465947[44] R. Zippel. 1979. Probabilistic algorithms for sparse polynomials. In
EUROSAMโ79(LNCS) , Vol. 72. Springer, 216โ226., Vol. 72. Springer, 216โ226.