<\body> <\doc-data|||<\author-affiliation> LIX, CNRS École polytechnique 91128 Palaiseau Cedex France ||> \; >|<\doc-note> This work has been supported by the ANR-09-JCJC-0098-01 project, as well as a Digiteo 2009-36HD grant and Région Ile-de-France. > \; >> In this paper, we study the complexity of several basic operations on linear differential operators with polynomial coefficients. As in the case of ordinary polynomials, we show that these complexities can be expressed almost linearly in terms of the cost of multiplication. ||> Let > be an effective field of constants of characteristic zero, so that all field operations can be carried out by algorithms. Given an indeterminate and the derivation =x*\>, where =\/\ x>, it is well known that the skew polynomial ring |]>> behaves very much like an ordinary polynomial ring: there are skew analogues for each of the classical operations of division with remainder, greatest common divisors, least common multiples, etc. In this paper, we will study the complexity of these operations. For this purpose, it will be more appropriate to work in the ring |]>> instead of |]>>. In analogy with the commutative case, we will give bounds for the computational complexities of the various operations in terms of the complexity of operator multiplication. For our complexity measures, we make the customary assumption that all field operations in can be carried out in constant time >. We will try to express the complexities of our algorithms in terms of the following standard complexities: <\itemize> The time > required for the multiplication of two polynomials of degrees>n> and coefficients in >. It is classical that =O> and =O*> if > admits sufficiently many >-th roots of unity. The complexity >|)>> of multiplying two r> matrices with entries in >. It is classical that \2.376>, although \3> in practice. We will denote by > the subset of > of polynomials of degree n>. Likewise, we denote by |]>> the set of operators \|]>> of degree L\n> in and degree > L\r> in >. Now consider two linear differential operators \|]>>. We start with studying the following complexities: <\itemize> The complexity > of multiplying and . The complexity > of applying to a vector of polynomials in >. The cost > to compute a fundamental system of solutions to the monic equation +L|)> f=0> in |]>>, up to order |)>>, while assuming the existence of such afundamental system. Given a vector of truncated power series in >, the cost >> of computing a monic operator in +\|]>> with =O|)>>. The special case was first studied in, where it was shown that =O>|)>>, using evaluation-interpolation techniques. The inverse bound >=O|)>> has been proved in; this paper also contains detailed information on the constant factors involved in these bounds. Recently (and after the writing of a first version of this paper), the quasi-optimal bound =-2>|)>> was proved in. For fixed constants ,\\0>, one has *n|)>=O|)>>, *r|)>>=O>|)>>, *n,\*r|)>>=O|)>>, etc., by splitting the multiplicands in a finite number of pieces. In this paper, we will freely use this remark without further mention. In order to simplify our complexity estimates, it will be convenient to make a few additional assumptions. First of all, we will assume that \2>, whence in particular *log n=O-1>|)>>. We will also assume that the function /n> is increasing and that />> is increasing in both and . This will indeed be the case for the complexity bounds for > that will be given in Section. In Section, we will first prove (see Theorems and) that the problems of multiplication and operator-vector application are essentially equivalent whenr>. We also recall the best existing bounds for operator multiplication. In Section, we show that the problems of computing fundamental systems of solutions and its inverse can be reduced to operator multiplication modulo a logarithmic overhead (see Theorems and). This provides a dual way to perform operations on differential operators by working on their fundamental systems of solutions. In Section and all subsequent sections, we always assume that r>. This is indeed required for the truncations of the fundamental systems of solutions at order |)>> to be linearly independent. In Section, we start with the operations of exact right division and right division with remainder. In Section, we consider greatest common right divisors (gcrds) and least common left multiples (lclms). Again, we will show how to express the complexities of these operations essentially in terms of the complexity >> of multiplication(see Theorems, , and ). For several of our algorithms, we need to work at a point where certain operators are non singular. If we only need the input operators to be non singular, then it is easy to find a point where this is the case. If we also need the output operators or certain auxiliary operators to be non singular (as in Section), then we resort to picking random points, which are non singular with probability1. In Section we present additional techniques for turning algorithms which rely on random point picking into randomized algorithms of Las Vegas type and into fully deterministic algorithms. For technical reasons, we found it convenient to work with respect to the Euler derivation > instead of >. Nevertheless, operators in |]>> can be converted efficiently into operators in |]>> and , modulo an increase of the degree in with the degree in > or > (see Lemma). Using our assumption that r>, such increases of the degree by only gives rise to constant overheads in the complexity bounds. Hence, the complexity bounds for our main algorithms from Sections, and still hold when replacing > by>. In addition, some of the algorithms can be adapted to directly use > instead of >, without the need for any conversions (see Remark). To the best of our knowledge, the idea to perform operations on linear differential operators power series solutions was first proposed (but only partially worked out) in. In this paper, we use a slightly different technique: instead of a single power series solution, we prefer to consider a fundamental system of solutions. This has the advantage of forcing a clean bijection between operators and solution spaces, thereby avoiding part of the randomness in the proposals from. It is also possible to mimic classical divide and conquer algorithms for right division, greatest common right divisors and least common left multiples, while using adjoints in order to perform the recursive operations on the appropriate side. Such algorithms were partially implemented inside and we plan to analyze this technique in more details in a forthcomingpaper. Various complexity results for computations with linear differential operators and other skew polynomials were previously obtained . Especially the computation of greatest common right divisors and least common left multiples of two or more operators has received particular attention. After the publication of a first version of this paper, the complexities of several classical algorithms for the computation of least common right multiples were studied in great detail in, and new improvements were proposed there. The complexities of most of the algorithms in this paper are stated in terms of the input output sizes. The uncertified randomized algorithms for gcrds and lclms are optimal up to logarithmic factors from this perspective, which yields an improvement with respect to the previously known complexity bounds. In the context of certified randomized algorithms ( of Las Vegas type), the complexity bounds remain quasi-optimal in terms of the size of a suitable certificate. From the deterministic point of view, the new algorithms for gcrds and lclms are suboptimal. The key argument behind the proof from that =O>|)>> is the observation that an operator \|]>> is uniquely determined by its images on the vector =,x|)>>>. This makes it possible to use a similar evaluation-interpolation strategy for the multiplication of differential operators as in the case of FFT-multiplication of commutative polynomials. More precisely, given \|]>>, let > be the matrix of the mapping \\;P\L> with respect to the bases > and >: <\eqnarray*> >||>|>||)>>>|>||>>|>|>||)>>>>>>.>>>> The evaluation and interpolation steps can be done efficiently using the following lemma, which is essentially contained in: <\lemma> Let \|]>>. Then <\enumerate-alpha> We may compute > as a function of in time *log r|)>>. We may recover from > in time *log r|)>>. <\proof> Consider the expansion of with respect to <\eqnarray*> |)>>|||)>+\+x*L|)>.>>>> For all , we have <\eqnarray*> |)>|)>>||*L|)>|]>|)>>>|||*L+j|)>|]>>>|||.>>>> In other words, > is a lower triangular band matrix <\eqnarray*> >||>||>|>|>|>|>||>>||>|>>|||>>>>>>>>> of bandwidth n>. The coefficients on the -th subdiagonal of > are exactly the result of amultipoint evaluation of > at ,r-1>. It is classical that both multipoint evaluation and the inverse operation of interpolation can be performed in time *log r|)>>. Doing this for each of the polynomials ,\,L> yields the result. <\theorem> If r>, then <\eqnarray*> >||+n**log r|)>>>>> <\proof> Let \|]>> and assume that we want to compute . We may evaluate |)>> in time ,2*r|)>>=O|)>>. We may also evaluate |)>|)>> in time >=O|)>>. Using Lemma, we may recover from |)>|)>> in time *log r|)>>. This completes the proof. <\theorem> If r>, then <\eqnarray*> >||+n**log r|)>.>>>> <\proof> Assume now that we are given |)>\\|]>>, as well as a vector ,\,V|)>>\\> and that we want to evaluate =|)>,\,K|)>|)>>. This is equivalent to evaluating the operator >=K-r|)>> at the vector *V>. It is classical that >> can be computed in time |)>>. Using Lemma, we may compute the unique operator \|]>> with |)>=x*V> in time **log r|)>=O*log r|)>>. We may next compute the product >*L> in time >=O|)>>. Lemma finally allows us to evaluate >*L> at > in time *log r|)>>, thereby yielding >. The above results immediately imply the bound =O>|)>> from by the computation of a product to the computation of a matrix product <\eqnarray*> >||*\.>>>> After the publication of a first version of this paper, the following quasi-optimal bound for> was established in3 and5>. <\theorem> <\enumerate-roman> For n>, we have =O-1>*r+n**log r|)>>. For r>, we have =O-1>*n+r**log n|)>>. The inverse bound >=O|)>> from can also be generalized: <\theorem> If r>, then the product of an n> matrix and an r> matrix with coefficients in > can be computed in time |)>>. <\proof> By the result from, the problem is equivalent to the computation of |n/r|\>> operators ,\,K> in |]>> with a fixed operator |]>>>. Setting +x*K+\+x>*K>, we may compute in time |)>>. We may directly read off the products *L,\,K*L> from the result. In this paper, we have chosen to work with respect to the derivation > instead of >. The following result from 3.3> can be used to efficiently convert between operators in |]>> and |]>> (in , we proved a somewhat weaker result which would also suffice for the purposes of this paper). We have written |]>> for the set of operators of degree n> in and degree r> in >. <\lemma> <\enumerate-alpha> Any operator in |]>> can be converted into an operator in |]>> in time **log r|)>>. Any operator in *\|]>> can be converted into an operator in |]>> in time **log r|)>>. From now on, we will assume that r>. We recall that an operator \|]>> of order is said to be at >, if its leading coefficient > does not vanish at >. We will say that an operator \|]>> of order is non singular (at the origin) if r>*L\|]>>> and *L> is non singular as an operatorin>. Given a non singular differential operator \|]>> of order , the equation =0> admits a fundamental system ,\,H|)>>> of solutions in |]>>, with the property that |)>=1> and |)>=0> for all r>> with j>. Conversely, given a >-linearly independent vector of power series \|]>>, there exists aunique monic operator \+\|]>|]>> of order with =0>. Let us show how to convert efficiently between these two representations. <\theorem> Let \|]>> be a differential operator of order n>, which is non singular at the origin, and let be its canonical fundamental system of solutions. Then we may compute up to order |)>> in time *log n|)>>. In other words, <\eqnarray*> >||*log n|)>.>>>> <\proof> Modulo multiplying on the left by 1>>, we may assume without loss of generality that is monic. Since is non singular at the origin, we have *L\\|]>>. Rewritten in terms of >, this means that is of the form <\eqnarray*> |||)>+x*C*\|)>+\+x*C*\|)>.>>||)>>||*-1|)>*\*-k+1|)>,>>>> for certain ,\,C\\>. Setting |)>-L\x*\|]>>, we observe that maps |]>> into *\|]>>. We now compute using the \Precursive\Q formula <\eqnarray*> ||>|>>|>>>>>+\|)>|)>,>>>> where <\eqnarray*> |)>r>A*x|)>>||r>|\>*x.>>>> The equation() is a schoolbook example for applying the strategy of relaxed resolution of power series equations. Since |)>> operates coefficientwise, it can be computed in linear time. The main cost of the computation therefore reduces to the relaxed evaluation of >. Using fast relaxed multiplication, this amounts to a cost <\eqnarray*> >|||n/2|\>,r|)>+4*|n/4|\>,r|)>+\+n*.>>>> Using the monotonicity assumption and Theorem, the result follows. In what follows, given a non zero series in , we denote by > its valuation. Given a vector of elements in a >-vector space, we will also denote by > the subvector space generated by the entries of , and <\eqnarray*> >||:Y\Vect\|}>.>>>> Notice that \dim|)>-1>. <\theorem> Let ,\,H|)>\\|]>> be >-linearly independent. Then there exists a unique monic operator \\+\|]>|]>> with =0>. Moreover, given the truncation of at order |)>>, we may compute at order >|)>> in time *log r|)>>. In other words, <\eqnarray*> >||*log r|)>.>>>> <\proof> Modulo a triangularization of , we may assume without loss of generality that |)>\\\v|)>=v>. We define operators >,\,L>> by <\eqnarray*> >>||>|>>||- L>|)>|L>|)>>|)>*L>.>>>> Then >> annihilates and for any other operator \\+\|]>|]>> with =0>, we would have -L|)>=0>, which is in contradiction with the fact that -L|)>\r>>. Moreover, by induction over , we observe that the coefficient of > in >> is given by -v|)>|)>*\*-v|)>|)>> and the coefficients of ,\,x> in >> can be expressed in terms of the coefficients of ,\,x|)>>> in ,\,H>. In particular, the truncation of at order >|)>> is uniquely determined by the truncation of at order |)>>. In order to explicitly compute up to a given order, it is more efficient to use a divide and conquer approach. More precisely, given ,\,H|)>\\> we compute \\+\|]>> using the following method: <\itemize> If , then we take =\- H/H|)> mod x>. Otherwise, let with |r/2|\>>. Compute ann,\,H|)>>. Evaluate |)>,\,A|)>|)> mod x>. Compute ann,\,I|)>>. Return >. If v>, then it is easy to check that =O>|)>>. For a fixed constant , we thus have <\eqnarray*> >|>|+C*.>>>> The result now follows from the monotonicity assumption. <\remark> If /r>> is increasing in for some \0>, then the bound further simplifies to =O|)>>. <\remark> We notice that the operator in Theorem is singular if and only if =r-1>, and if and only if :Y\Vect\|}>=,r-1|}>>. <\remark> The algorithm from the proof can be adapted so as produce a vanishing operator in *\+\|]>|]>> instead of +\|]>|]>>. Indeed, for this, it suffices to take <\eqnarray*> >>||- L>|)>|L>|)>>|)>*L>,>>>> and carefully adapt the truncation orders. Although a general operator \|]>> can be singular at the origin, many operations on operators (such as right division and greatest common right divisors) commute with translations x+x>, and Lemma may be used in conjunction with the following lemma in order to reduce to the case when is non singular at the origin. <\lemma> Given a non zero operator \|]>>, we may find a point \\> where is non singular in time |)>>. <\proof> Let > be the leading coefficient of . Since L\n>, we have |)>\0> for some \,n|}>>. Using fast multipoint evaluation, we may find such a point > in time |)>>. From the formula() it is clear that both the degrees in and > are additive for the multiplication of operators \|]>>>. In particular, if \|]>> and is left or right divisible by , then the quotient is again in |]>>. <\theorem> Let \|]>> be such that for some \|]>> and r>. Then we may compute in time *log n|)>>. <\proof> By Lemmas and, and modulo ashift x+x>, we may assume without loss of generality that and are non singular at the origin. We now use the following algorithm: <\itemize> We first compute the canonical fundamental system of solutions to =0>> up to order |)>>. By Theorem , this can be done in time *log n|)>>>. We next evaluate > and compute a>-basis for > at order |)>>. This can be done in time |)>>, by Theorems and, and using linear algebra. Since is non singular, we have \deg> K\v|)>=v> for all \|]>>. In particular, the > Q=deg> L-deg> K> elements of of valuations > K,\,deg> L-1> are mapped to set which spans a vector space of dimension > Q>. This shows that mod x|)>=deg> Q>. We now compute the monic annihilator =ann> of at order |)>>. This can be done in time *log r|)>=O*log n|)>>, by Theorem. We return the truncation of *\> at order |)>>, where =L> L>/K> K>>. Since each of the steps can be carried out in time *log n|)>>, the result follows. It is classical that euclidean division generalizes to the skew polynomial ring |]>>. In other words, given operators \|]>> where 0>, there exist unique operators > and > in |]>> with <\eqnarray*> ||>>> and > R\deg> B>. If \|]>> and is the leading term of with respect to >, then left multiplication of by > A-deg> B+1>> allows us to remain in the domain |]>>>: there exist unique > and > in |]>> with <\eqnarray*> > A-deg> B+1>*A>||>>>> and > R\deg> B>. The operators and are usually called pseudo-quotients and pseudo-remainders. In some cases, a non trivial polynomial can be factored out in the relation(). Let be monic, of maximal degree, such that *Q*B,J1>*R\\|]>>. Then we call 1>*Q=quo>> and 1>*R=rem>> the \Psimplified\Q pseudo-quotient and pseudo-remainder of and . <\lemma> Let ,\,H|)>\\|]>> be >-linearly independent and define |)>+1>. Given *\|]>|)>>, there exists a unique operator \|]>|]>> of order r> with =G> and we may compute its first terms with respect to in time *log n|)>>. <\proof> Let =v|)>> for each . Modulo a base change, we may assume without loss of generality that \\\\>. Let :\|]>\\|]>> be the operator with <\eqnarray*> ,\,V|)>>||>*V,\,x>*V|)>,>>>> and let > denote the inverse operator. Let :\|]>|]>\\|]>> be the operator with <\eqnarray*> >||1>|)>|)>.>>>> Writing K*x*\> and =|)>>, we have <\eqnarray*> >>|>>|>>>>>>|||>|>||)>>>|>|>||>>||>|>||)>>>>>>>*>>|>>|>>>>>.>>>> In other words, > and its inverse > operate coefficientwise and coefficients can be computed in time >*n|)>>. Putting =x>+E> with =o>|)>> for each , we may rewrite the equation =G>as <\eqnarray*> ||1>|)>|)>>>>> and we observe that the coefficient of > in the righthand side of() only depends on earlier coefficients of ,x> in . In particular, we may solve the equation using a relaxed algorithm. Then the main cost is concentrated in the relaxed evaluation of >. As in the proof of Theorem, this evaluation can be done in time *log n|)>>. <\theorem> Let \|]>> with r> and > K\0>. Right pseudo-division of by and simplification yields a relation <\eqnarray*> ||>>> where >,R=rem>\\|]>>. If \n> is such that \|]>,r>>>, then and can be computed in time ,r|)>*log n|)>>. <\proof> Modulo a shift x+x>, we may assume without loss of generality that and are non singular at the origin. We now use the following algorithm: <\itemize> We compute the canonical fundamental system of solutions to =0> up to order+r>|)>>. This requires a time ,s|)>*log n|)>>. We compute > with =A*G> up to order +r>|)>>. This requires atime ,r|)>|)>>. We determine the operator \\|]>|]>> with =x*G> up to order +r>|)>>. Lemma shows that this can be done in time ,s|)>*log n|)>>. By Theorem, we have *A*\> and s>*\> is known up to order >|)>>. Now s>*\,\,x*\> are truncated rational functions, for which the degrees of the numerators and denominators are bounded by >. Using rational function reconstruction, we may thus compute /D=x*\> with ,D|)>=1> in time *log n|)>>>. Taking ,\,D|)>>, we find . Once and are known, we compute using the algorithm from Theorem. The total complexity of this algorithm is bounded by ,r|)>*log n|)>>. <\remark> In the above proof, we have assumed that > is known beforehand. In general, we may still apply the above algorithm for a trial value>>. Then the algorithm may either fail (for instance, if ,\,D|)>\n>>), or return the triple > under the assumption that \|]>>,r>>. We may then check whether the triple is correct in time >,r|)>|)>>. Applying this procedure for successive guesses >=n,2*n,4*n,\>, the algorithm ultimately succeeds for an >> with >\2*n>. Using the monotonicity hypothesis, the total running time thus remains bounded by >,r|)>*log n>|)>=O,r|)>*log n|)>>. It is classical that greatest common right divisors and least common left multiples exist in the skew euclidean domain |]>>: given two operators \|]>>>, the greatest common right divisor =gcrd> and the least common left multiple =lclm> are the unique monic operators with <\eqnarray*> |]>*\>|||]>*K+\|]>*L>>||]>*\>|||]>*K\\|]>*L.>>>> Assume now that \|]>> and let and be monic polynomials of minimal degrees, such that > and > are in |]>>. Then we call >=gcrd>=A*\> and >=>>=B*\> the (simplified) pseudo-gcrd and pseudo-lclm of and . <\lemma> Let \|]>> be such that and are non singular at the origin, as well as >> or >>. Let and be the canonical fundamental systems of solutions to =0> and =0>. Then <\eqnarray*> > lclm>>||+Vect|]> mod x|)>>>|> gcrd>>||\Vect|]> mod x|)>.>>>> <\proof> Let >=gcrd>>, >=lclm>>, > \>> and > \>\2*r>. If >> is non singular, then it admits a canonical fundamental system of solutions ,\,M|)>> with |)>=1> and |)>=0> for all t> with j>. In particular, mod x|)>=t>. Since >> is the least common left multiple of and , we also have =Vect+Vect>, which completes the proof of the first equality. If >> is non singular, then we obtain the second equality in a similar way. If >> is non singular, then we also have mod x|)>=deg> K> and mod x|)>=deg> L>, since and are non singular. Now \Vect|]> mod x|)>=dim mod x|)>+dim mod x|)>-dim+Vect|]> mod x|)>>, whence \Vect|]> mod x|)>=deg> K+deg> L-deg> \>=deg> \>>. If >> is non singular, then we obtain the first equality in a similar way. <\theorem> Let \|]>> and \n> be such that >=gcrd>\\|]>,r>> and r>. Assume that , and >> (or >>) are non singular at the origin. Then we may compute>> in time ,r|)>*log n|)>>. <\proof> We compute >> using the following algorithm: <\itemize> We compute the canonical fundamental systems of solutions and to =0> and =0> at order +r>|)>>. This can be done in time ,r|)>*log n|)>>. Using linear algebra, we compute a basis for \Vect> at order +r>|)>>. This can be done in time n*r-1>|)>>. By Lemma, we have dim|)>=deg> \>>. We also notice that \r>. We compute =ann=gcrd> at order >|)>>. By Theorem, this can be done in time ,r|)>*log n|)>>. We compute >> from mod x>> using rational function reconstruction. This algorithm requires a total running time ,r|)>*log n|)>>. <\remark> In the above proof, we have again assumed that > is known beforehand. Below, we will discuss ways to check the correctness of the computed result for a trial value>>, after which a similar strategy as in remark can be applied. During the relaxed computation of and , we may also check whether > at each next coefficient. In the particular case when =1>, the running time of the algorithm will then be bounded by >,r|)>*log n>|)>>, where >> is the smallest order at which common solutions no longer exist. This kind of early termination only works for this very special case. <\remark> Notice that >> might be singular at the origin, even if , and >>> are not. This happens for instance when is the minimal annihilator of the vector > and the minimal annihilator of the vector ,x|)>>, so that =\-1>. <\theorem> Let \|]>> and \n> be such that >=lclm>\|]>,2*r>>> and r>. If , and >> (or >>) are non singular at the origin, then we may compute>> in time ,r|)>*log n|)>>. <\proof> Similar to the proof of Theorem, by taking +Vect> instead of \Vect>. The assumption that >> should be non singular is still a bit unsatisfactory in Theorems and, even though the probability that a randomly chosen point is singular is infinitesimal. \ If we drop this assumption, then we still have deg> \>> in the proof of Theorem. Consequently, \Pcandidate\Q pseudo-gcrds >> found by the algorithm are genuine pseudo-gcrds whenever >> pseudo-divides both and . Using the right division algorithms from the previous section, this can be checked in time *r,r|)>*log n|)>> in the case of gcrds and *log n|)>> in the case of lclms. <\remark> Using the polynomial linear algebra techniques from, it is likely that one may prove that >> for some \> and \|]>>. If this is indeed the case, then the trial divisions of and by >> can actually be carried out in time *log n|)>>. An alternative way to check whether candidate gcrds and lclms are correct is to compute Bezout and Ore relations. More precisely, given \|]>> with \*K>, there exist operators \\|]>> with <\eqnarray*> >>|>>>>>|||>||>>>>*>|>>>>,>>>> > A*K,deg> B*L\deg> \> and >. The 2> matrix at the righthand side will be called the Euclidean matrix =Eucl> of and . In a similar way as above, we may define a (simplified)pseudo-Euclidean matrix >=Eucl>> with entries >,B>,C>,D>>> in |]>>, whenever \|]>>. We will say that > is non singular at >, if the denominators of and do not vanish at >. <\theorem> Let \|]>> and \n> be such that >=Eucl>\|]>,r>2>>> and r>. If , , >> and > are non singular at the origin, then we may compute >> in time ,r|)>*log n|)>=*r-1>|)>>. <\proof> Assuming > known, we compute =|>||>>>>> at order >|)>> as follows: <\itemize> We compute the canonical fundamental systems of solutions and to =0> and =0> at order +3*r>|)>>. We compute a basis for \Vect> at order +3*r>|)>>, together with bases > and > for the supplements of > in > >. We also compute =ann> at order +2*r>|)>>. We solve the systems |)>|)>=\|)>> and |)>|)>=\|)>> in at order >|)>>, using Lemma. We compute a basis for +Vect> at order +2*r>|)>>, as well as bases > and > for the vector spaces |)>> |)>> at order +2*r>|)>>. We compute > K>*ann|)>> and > L>*ann|)>> at order >|)>>. We finally compute >> from and using rational function reconstruction. The complexity analysis and the remainder of the proof is done in a similar way as in the proofs of Theorems and. With the above techniques, we may at least verify whether computed pseudo-gcrds or pseudo-lclms are correct. For a fully deterministic algorithm, we still need away to find a point where >> is non singular. This can be done by brute force. Let us state the result in the most general setting of pseudo-Euclidean matrices. <\theorem> Let \|]>> and \n> be such that >=Eucl>\|]>,r>2>>> and r>. Then we may compute >> in time ,r|)>*log n+n**r+r>*log r|)>|)>=*>+n*r|)>|)>>. <\proof> Let > K>, > L>, and assume first that we know >. Let ,\,x+n>> be +n+1> be pairwise distinct, randomly chosen points in > at which and are non singular. At each >, we compute canonical fundamental systems of solutions and for and at order |)>>. We claim that this can be done in time |)>*r*log n+n**r+r>*log r|)>|)>>. Indeed, it requires a time **log r|)>> to rewrite each operator with respect to>. We next perform a multipoint evaluation of the coefficients of these operators to obtain the shifted operators at ,\,x+n>> (this requires a time |)>*r*log n|)>>). The truncations of these operators at order |)>> are then converted back to the respresentation with respect to >. This can be done in time *r**log r|)>>. Using Theorem, we finally compute the required fundamental systems of solutions in time **log r|)>=O*r>*log r|)>>. From >\\|]>,r>2>>, we get >=lclm>\\|]>+n,2*r>>. Since we assumed > to be sufficiently large, it follows that >=lclm>> is non singular at one of the points>. At such a point >, the canonical fundamental systems of solutions and generate avector space +Vect> of maximal dimension deg> \>>, and with a basis ,\,y> such that |)>=k> for all k\s>. We finally apply Theorem in order to obtain >>. If > is unknown, then we use a sequence of guesses =n,2*n,4*n,\>, as in the previous proofs. <\remark> In the case of least common left multiples, we may directly compute >> using Theorem and certify the result using trial division by and . This allows us to use the weaker assumption >\\|]>,2*r>> instead of >\\|]>,r>2>>, whereas the complexity bound becomes *log n+n**r+r>*log r|)>|)>=*>+n*r|)>|)>>. We have summarized our complexity bounds for Euclidean operations on two operators \|]>> in Table. We systematically write > for the degree in of the result. We also write >> for the degree of the Euclidean matrix in . The algorithms in the first line correspond to applying Theorems, and at arandomly chosen point, without checking the result. The second line corresponds to the Las Vegas randomized algorithm for which the answers are certified through trial division (the bound for gcrds might further drop to >|)>> in view of Remark; more generally, the bounds can be restated in terms of sizes of certificates). In the third line, we rather use Euclidean matrices for the certification. The fourth line shows complexity bounds for the deterministic versions of our algorithms. In comparison, several randomized Las Vegas algorithms were given in that achieve the complexity bound >|)>> for lclms. This is in particular the case for Heffter's algorithm, when using Theorem. The non determinism is due to the use of a fast Las Vegas randomized algorithm for the computation of kernels of matrices with polynomial entries2>. Grigoriev established complexity bounds for gcrds which rely on asimilar reduction to polynomial linear algebra. Along the same lines as in, this should lead to a Las Vegas randomized algorithm of complexity >|)>>, although we did not check this in detail. In summary, the new algorithms do not achieve any improvements in the worst case. Nevertheless, the uncertified versions of our algorithms admit optimal running times up to logarithmic factors in terms of the combined input output size. The certified randomized versions satisfy similar complexity bounds in terms of the size of a suitable certificate; such bounds can sometimes be better than the previously known worst case bounds. When performing our expansions at arandomly chosen point in >, we also recall that the probability of failure is exponentially small as a function of the bitsize of this point. <\big-table|||||>||*r-1>|)>>>|*r-1>|)>>>|*r-1>|)>>>>| division>|*r>|)>>>|>|)>>>|>||>*r-1>|)>> >|>*r-1>|)>> >|*r-1>|)>>>>||>*>+n*r|)>|)>>>|*>+n*r|)>|)>>>|*>+n*r|)>|)>>>>>>>> Complexity bounds for the Euclidean operations on two operators and . The algorithms from Section extend in a straightforward way to the computation of greatest common right divisors and least common left multiples of more than two operators. For instance, using obvious notations, we obtain the following generalizations of Theorems and. <\theorem> Let ,\,L\\|]>> with r> and \r>, \max|)>> be such that >=>,\,L|)>>\\|]>,r>>>. Assume that ,\,L> and >> are all non singular at the origin. Then we may compute>> in time ,r|)>*log n+,r|)>*log n>+k*r*|)>-2>*n|)>>. <\proof> We compute >> using the following algorithm: <\itemize> We first compute the canonical fundamental systems of solutions > to |)>=0> at order +r>|)>>. This can be done in time ,r|)>*log n|)>>. Let =Vect|)>+\+Vect|)>> for all i\j\k>. Using linear algebra, we may recursively compute a basis > for > from bases > and > for > and >, where |/2|\>>. This algorithm yields a basis for > in time |)>-2>*n|)>>. Using a suitable generalization of Lemma, we also notice that >|)>=deg> \>>. We compute =ann=lclm,\,L|)>> at order >|)>>. By Theorem, this can be done in time ,r|)>*log n|)>>. We compute >> from mod x>> using rational function reconstruction. We obtain the result by adding up all complexity bounds. <\remark> When taking =k*r\n> and using, the complexity bound simplifies to ,k*r|)>*log n|)>=O-1>*r-1>*n*log n+k*r*|)>*log n|)>>. By 6>, we may always take =n*r*k>, after which the bound further reduces to +1>*r>*n|)>>. In our randomized setting, this improves upon the bounds from1>. <\remark> \ If we also require a certification of the result, then we may use the trial division technique. This amounts to exact divisions of operators in |]>+n*r,r>>> by ,\,L>>. Using the division algorithm from Section, and taking =k*r\n> and =n*r*k> as above, this can be done in time <\equation*> O+n*r,r|)>*log |)>|)>=+n*r|)>*|)>-1>|)>=+2>*r>*n|)>. This is slightly better than the new bound from. <\theorem> Let ,\,L\\|]>> and \n\r> be such that >=>,\,L|)>>\\|]>,r>>>. Assume that ,\,L> and >> are all non singular at the origin. Then we may compute>> in time ,r|)>*log n>+*log r>|)>>. <\proof> The proof is similar to the one of Theorem, except for the way how we compute a basis for |)>\\\Vect|)>>. Indeed, we first compute a basis > for >. This requires a time *log r|)>> for the computation of ,\,H> modulo > and a time >|)>> for the remaining linear algebra. We next compute the unique constant matrix such that > modulo >. Since >> is non singular, we have > at any order, so it suffices to compute > up to order +r>> in order to obtain up to order +r>>. <\remark> An interesting question is whether there exists a faster algorithm to compute the orders and of >=gcrd>,\,L|)>> and >=lclm>,\,L|)>>, without computing >> and >> themselves. For this, it suffices to compute the dimensions of |)>\\\Vect|)>> and |)>+\+Vect|)>>. Assuming that we are at a \Pnon singular point\Q, the answer is therefore yes: using the techniques from the proofs of Theorems and, we may compute in time *log r|)>=>|)>> and in time *log t+k*r*t-1>|)>>=-1>|)>>. <\bibliography|bib|tm-plain|all.bib> <\bib-list|35> A.V.Aho, K.SteiglitzJ.D.Ullman. Evaluating polynomials on a fixed set of points. , 4:533\U539, 1975. A.Benoit, A.BostanJ.vander Hoeven. Quasi-optimal multiplication of linear differential operators. , 524\U530. New Brunswick, October 2012. IEEE. A.BorodinR.T.Moenck. Fast modular transforms. , 8:366\U386, 1974. A.Bostan. . , École polytechnique, 2003. A.Bostan, F.ChyzakN.Le Roux. Products of ordinary differential operators by evaluation and interpolation. J.Rafael SendraLaureano González-Vega, , 23\U30. Linz/Hagenberg, Austria, July 2008. ACM Press. A.Bostan, F.Chyzak, B.SalvyZ.Li. Fast computation of common left multiples of linear ordinary differential operators. J.vander HoevenM.van Hoeij, , 99\U106. Grenoble, France, July 2012. A.BostanÉ.Schost. Polynomial evaluation and interpolation on special sets of points. , 21(4):420\U446, August 2005. Festschrift for the 70th Birthday of Arnold Schönhage. E.Brassine. , 331\U347. Note III du Tome 2 du Cours d'analyse de Ch. Sturm. École polytechnique, 1864. D.G.CantorE.Kaltofen. On fast multiplication of polynomials over arbitrary algebras. , 28:693\U701, 1991. J.W.CooleyJ.W.Tukey. An algorithm for the machine calculation of complex Fourier series. , 19:297\U301, 1965. D.CoppersmithS.Winograd. Matrix multiplication via arithmetic progressions. > Annual Symposium on Theory of Computing>, 1\U6. New York City, may 25\U27 1987. S.S.Demidov. On the history of the theory of linear differential equations. , 28(4):369\U387, 1983. F.Le Gall. Powers of tensors and fast matrix multiplication. , 296\U303. Kobe, Japan, July 23\U25 2014. M.Giesbrecht. Factoring in skew polynomial rings over finite fields. , 583, 191\U203. 1992. M.Giesbrecht. Factoring in skew polynomial rings over finite finite fields. , 26:463\U486, 1998. M.GiesbrechtY.Zhang. Factoring and decomposing ore polynomials over >. Manuel Bronstein, , 127\U134. Philadelphia, USA, August 2003. D.Yu.Grigoriev. Complexity of factoring and calculating the gcd of linear ordinary differential operators. , 10(1):7\U37, 1990. J.vonzur GathenJ.Gerhard. . Cambridge University Press, 2-nd, 2002. L.Heffter. Über gemeinsame Vielfache linearer Differentialausdrücke und lineare Differentialgleichungen derselben Klasse. , 116:157\U166, 1896. J.vander Hoeven. FFT-like multiplication of linear differential operators. , 33(1):123\U127, 2002. J.vander Hoeven. Relax, but don't be too lazy. , 34:479\U542, 2002. J.vander Hoeven. Relaxed multiplication using the middle product. Manuel Bronstein, , 143\U147. Philadelphia, USA, August 2003. J.vander Hoeven. On the complexity of skew arithmetic. , HAL, 2011. . J.vander Hoeven, G.Lecerf, B.Mourrain etal. Mathemagix. 2002. . Z.Li. A subresultant theory for ore polynomials with applications. O.Gloor, , 132\U139. Rostock, Germany, August 1998. G.Libri. Mémoire sur la résolution des équations algébriques dont les racines ont entre elles un rapport donné, et sur l'intégration des équations différentielles linéaires dont les intégrales particulières peuvent s'exprimer les unes par les autres. , 10:167\U194, 1833. R.T.MoenckA.Borodin. Fast modular transforms via division. , 90\U96. Univ. Maryland, College Park, Md., 1972. O.Ore. Theorie der linearen Differentialgleichungen. , 167:221\U234, 1932. O.Ore. Theory of non-commutative polynomials. , 34(3):480\U508, 1933. V.Pan. , 179 Springer, 1984. A.Schönhage. Schnelle Multiplikation von Polynomen über Körpern der Charakteristik 2. , 7:395\U398, 1977. A.SchönhageV.Strassen. Schnelle Multiplikation groÿer Zahlen. , 7:281\U292, 1971. R.P.Stanley. Differentially finite power series. , 1:175\U188, 1980. MR #81m:05012. V.Strassen. Gaussian elimination is not optimal. , 13:352\U356, 1969. V.Strassen. Die Berechnungskomplexität von elementarsymmetrischen Funktionen und von Interpolationskoeffizienten. , 20:238\U251, 1973. <\initial> <\collection> <\attachments> <\collection> <\associate|bib-bibliography> <\db-entry|+2UJhhqqcr6IfEs|article|Lib1833> <|db-entry> > <\db-entry|+2UJhhqqcr6If9G|inbook|Bras1864> <|db-entry> > <\db-entry|+2UJhhqqcr6IfAT|article|Dem83> <|db-entry> > <\db-entry|+2UJhhqqcr6IfGc|article|Ore32> <|db-entry> > <\db-entry|+2UJhhqqcr6IfGd|article|Ore33> <|db-entry> > <\db-entry|+2UJhhqqcr6IfJ7|article|SS71> <|db-entry> V. > <\db-entry|+2UJhhqqcr6IfIF|article|Sch77> <|db-entry> > <\db-entry|+2UJhhqqcr6If9w|article|CK91> <|db-entry> E. > <\db-entry|+2UJhhqqcr6IfAB|article|CT65> <|db-entry> J.W. > <\db-entry|+2UJhhqqcr6IfJJ|article|Str69> <|db-entry> > <\db-entry|+2UJhhqqcr6IfGl|book|Pan84> <|db-entry> > <\db-entry|+2UJhhqqcr6IfAC|inproceedings|CW90> <|db-entry> S. > > Annual Symposium on Theory of Computing> <\db-entry|+2UJhhqqcr6IfBm|inproceedings|Gall14> <|db-entry> > <\db-entry|+2UJhhqqcr6IfLG|article|vdH:FFT> <|db-entry> > <\db-entry|+2UJhhqqcr6If8B|inproceedings|BCR08> <|db-entry> F. N. Le > Laureano > <\db-entry|+2UJhhqqcr6IfMc|inproceedings|vdH:ldomul> <|db-entry> A. J. van der > <\db-entry|+2UJhhqqcr6If93|phdthesis|Bos:phd> <|db-entry> > <\db-entry|+2UJhhqqcr6IfN7|misc|vdH:mmx> <|db-entry> G. B. > > <\db-entry|+2UJhhqqcr6IfCT|article|Gri90> <|db-entry> > <\db-entry|+2UJhhqqcr6IfBy|inproceedings|Gies92> <|db-entry> > <\db-entry|+2UJhhqqcr6IfBz|article|Gies98> <|db-entry> > <\db-entry|+2UJhhqqcr6IfEr|inproceedings|Li98> <|db-entry> > > <\db-entry|+2UJhhqqcr6IfC0|inproceedings|GiZh03> <|db-entry> Y. > >> > <\db-entry|+2UJhhqqcr6IfML|techreport|vdH:skewcompl:pre> <|db-entry> > > <\db-entry|+2UJhhqqcr6IfCy|article|Heff1896> <|db-entry> > <\db-entry|+2UJhhqqcr6IfJC|article|Stan80> <|db-entry> > <\db-entry|+2UJhhqqcr6If8D|inproceedings|BCSL12> <|db-entry> F. B. Z. > M. van > <\db-entry|+2UJhhqqcr6IfFR|inproceedings|MB72> <|db-entry> A. > <\db-entry|+2UJhhqqcr6IfJK|article|Str73> <|db-entry> > <\db-entry|+2UJhhqqcr6If8q|article|BM74> <|db-entry> R.T. > <\db-entry|+2UJhhqqcr6If7j|article|ASU75> <|db-entry> K. J. D. > <\db-entry|+2UJhhqqcr6IfL1|article|vdH:relax> <|db-entry> > <\db-entry|+2UJhhqqcr6IfLK|inproceedings|vdH:issac03> <|db-entry> > > <\db-entry|+2UJhhqqcr6If94|article|BoSc05> <|db-entry> É. > <\db-entry|+2UJhhqqcr6IfBk|book|GaGe02> <|db-entry> J. > <\references> <\collection> > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > <\auxiliary> <\collection> <\associate|bib> Lib1833 Bras1864 Dem83 Ore32 Ore33 SS71 Sch77 CK91 CT65 Str69 Pan84 CW90 Gall14 vdH:FFT BCR08 vdH:ldomul Bos:phd Bos:phd vdH:mmx Gri90 Gies92 Gies98 Li98 GiZh03 Bos:phd vdH:skewcompl:pre Heff1896 Stan80 Li98 BCSL12 vdH:FFT BCR08 MB72 Str73 BM74 ASU75 vdH:FFT vdH:ldomul BCR08 BCR08 BCR08 vdH:FFT vdH:relax vdH:issac03 BoSc05 GaGe02 Gri90 BCSL12 BCSL12 Heff1896 BCSL12 BCSL12 vdH:ldomul BCSL12 BCSL12 BCSL12 <\associate|table> <\tuple|normal> Complexity bounds for the Euclidean operations on two operators |K> and |L>. > <\associate|toc> |math-font-series||font-shape||1.Introduction> |.>>>>|> |math-font-series||font-shape||2.Evaluation and interpolation> |.>>>>|> |math-font-series||font-shape||3.Local solutions> |.>>>>|> |math-font-series||font-shape||4.Right division> |.>>>>|> |math-font-series||font-shape||5.Euclidean operations> |.>>>>|> |5.1.Randomized algorithms |.>>>>|> > |5.2.Certifying correctness |.>>>>|> > |5.3.Summary of the complexity bounds for Euclidean operations |.>>>>|> > |5.4.Generalizations |.>>>>|> > |math-font-series||font-shape||Bibliography> |.>>>>|>