A Bostan
email: alin.bostan@inria.fr

T Krick
email: krick@dm.uba.ar

A Szanto
email: aszanto@ncsu.edu

M Valdettaro

Subresultants of (x -α) m and (x -β) n , Jacobi polynomials and complexity

Keywords: Subresultants, algorithms, complexity, Jacobi polynomials. 2010 13P15, 15B05, 33C05, 33C45, 33F10, 68W30

In an earlier article together with Carlos D'Andrea [BDKSV2017], we described explicit expressions for the coefficients of the order-d polynomial subresultant of (x -α) m and (x -β) n with respect to Bernstein's set of polynomials {(x -α) j (x -β) d-j , 0 ≤ j ≤ d}, for 0 ≤ d < min{m, n}. The current paper further develops the study of these structured polynomials and shows that the coefficients of the subresultants of (x -α) m and (x -β) n with respect to the monomial basis can be computed in linear arithmetic complexity, which is faster than for arbitrary polynomials. The result is obtained as a consequence of the amazing though seemingly unnoticed fact that these sub-

Introduction

Let K be a field, and let f = f m x m + • • • + f 0 and g = g n x n + • • • + g 0 be two polynomials in K[x] with f m = 0 and g n = 0. Set 0 ≤ d < min{m, n}. The order-d subresultant Sres d (f, g) is the polynomial in K[x] defined as Sres d (f, g) := det m+n-2d

f m • • • • • • f d+1-(n-d-1) x n-d-1 f . . . . . . . . . n-d f m . . . f d+1 f g n • • • • • • g d+1-(m-d-1) x m-d-1 g . . . . . . . . . m-d g n • • • g d+1 g , (1)
where, by convention, f = g = 0 for < 0. The polynomial Sres d (f, g) has degree at most d, and each of its coefficients is equal to a minor of the Sylvester matrix of f and g. In particular the coefficient of x d , called the principal subresultant of f and g, is given by PSres d (f, g) := det m+n-2d

f m • • • • • • f d-(n-d-1) . . . . . . n-d f m • • • f d g n • • • • • • g d-(m-d-1) . . . . . . m-d g n • • • g d .
Subresultants were introduced implicitly by Jacobi [Jac1836] and explicitly by Sylvester [Syl1839, Syl1840]; we refer to [Loo1983] and [GL2003] for detailed historical accounts 1 Let M(n) denote the arithmetic complexity of degree-n polynomial multiplication in K[x]. Precisely, M(n) is an upper bound for the total number of additions/subtractions and products/divisions in the base field K that are sufficient to compute the product of any two polynomials in K[x] of degree at most n. It is classical, see e.g. [GG2013, Ch. 8], that M(n) = O(n log n log log n) by using FFT-based algorithms. For arbitrary polynomials f, g ∈ K[x] of degree n, the fastest known algorithms are able to compute in O(M(n) log n) arithmetic operations in K either one selected polynomial subresultant Sres d (f, g) [Rei1997, LR2001, Lec2018], or all their principal subresultants PSres d (f, g) for 0 ≤ d < n [GG2013, Cor. 11.18]. It is an open question whether this can be improved to O(M(n)), even for the classical resultant (the case d = 0).

In this paper we present algorithms with linear complexity for these two tasks for the special family of polynomials considered in [BDKSV2017], namely f = (x -α) m and g = (x -β) n in K[x], when char(K) = 0 or char(K) ≥ max{m, n}, and α = β ∈ K (note that when α = β there is nothing to compute since all subresultants vanish). To our knowledge, we are exhibiting the first family of "structured polynomials" for which subresultants (and all principal subresultants) can be computed in optimal arithmetic complexity.

Let us first observe that the resultant Sres 0 ((x -α) m , (x -β) n ) = (αβ) mn , which corresponds to the case d = 0, can be computed by binary powering in O(log(mn)) arithmetic operations in K. The general case is not so simple: for example the particular case d = 1 of [BDKSV2017, Theorems 1.1 1 The Sylvester matrix was defined in [Syl1840], and the order-d subresultant was introduced in [Syl1839, Syl1840] under the name of "prime derivative of the d-degree". The term "polynomial subresultant" was seemingly coined by Collins [Col1967], and probably inspired to him by Bôcher's textbook [Boc1907, §69] who had used the word "subresultants" to refer to determinants of certain submatrices of the Sylvester matrix. Almost simultaneously, Householder and Stewart [HS1967, Hou1968] employed the term "polynomial bigradients". The principal subresultants were named "Nebenresultanten" (minor resultants) by Habicht [Hab1948]. The current terminology principal subresultants seems to appear for the first time in Collins' paper [Col1974]. and 1.2] (see also Theorem 2 below) shows that, for 1 < min{m, n},

Sres 1 ((x -α) m , (x -β) n ) = (α -β) (m-1)(n-1) m + n -2 m -1 x - m + n -3 m -1 α - m + n -3 n -1 β .
This identity implies that, from a computational perspective, there is already a striking difference between the cases d = 0 and d = 1. Indeed, although the term (α -β) (m-1)(n-1) can be computed in O(log(mn)) operations in K, no algorithm with arithmetic complexity polynomial in log(mn) is known for computing binomial coefficients such as m+n-2 m-1 . However, the right-hand side of the previous identity can be computed in O(min{m, n}) operations (see Lemma 8 below), provided the characteristic of the base field K is zero or large enough. The main result of the current article extends this complexity observation to arbitrary 1 ≤ d < min{m, n}.

Theorem 1. Let d, m, n ∈ N with 1 ≤ d < min{m, n} and let K be a field with char(K) = 0 or char(K) ≥ max{m, n}, and α, β

∈ K with α = β. Set Sres d ((x -α) m , (x -β) n ) = d k=0 s k x k .
Then, (a) if char(K) = 0 or char(K) ≥ m + n -d, then s d = 0 and all the coefficients s k for 0 ≤ k ≤ d can be computed using O(min{m, n} + log(mn)) arithmetic operations in K, (b) when char(K) = m + n -d -1, the following equality holds in K:

Sres d ((x -α) m , (x -β) n ) = (-1) md (α -β) (m-d)(n-d)+d and Sres d ((x -α) m , (x -β) n ) can be computed using O(log(mn)) arith- metic operations in K, (c) if m + n -d -1 > char(K) ≥ max{m, n}, then Sres d ((x -α) m , (x -β) n ) = 0.
We prove Theorem 1 via an amazing (and seemingly previously unobserved) close connection of the subresultants Sres d ((x-α) m , (x-β) n ) with the classical family of orthogonal polynomials known as the Jacobi polynomials, introduced and studied by Jacobi in his posthumous article [Jac1859]. This allows us to produce a recurrence for the coefficients of the subresultant, which is derived from the differential equation satisfied by the Jacobi polynomial, and hence by the subresultant.

To express the polynomial subresultants Sres d ((x-α) m , (x-β) n ) as Jacobi polynomials, let us recall [Sze1975, Chapter 4] that for any k, , r ∈ Z with r ≥ 0, the Jacobi polynomial

P (k, ) r (x) can be defined in 1 2 Z[x]
, and thus also in K[x] for any abstract field K with char(K) = 2, in two equivalent ways:

• by Rodrigues' formula

P (k, ) r (x) := (-1) r 2 r r! (1 -x) -k (1 + x) -∂ r ∂x r (1 -x) k+r (1 + x) +r ,
• as a hypergeometric sum:

P (k, ) r (x) := r j=0 (k + r -j + 1) j j! ( + j + 1) r-j (r -j)! x -1 2 r-j x + 1 2 j ,
where for any a ∈ Z, (a) 0 := 1 and (a) j := a(a + 1) • • • (a + j -1) for j ≥ 1 denotes the jth Pochhammer symbol, or the rising factorial, of a.

Our next result asserts that the d-th subresultant of (x -α) m and (x -β) n coincides, up to an explicit multiplicative constant and up to an affine change of variables, with the Jacobi polynomial P (-n,-m) d (x). More precisely, for α = β, we consider the following change of variables in the Jacobi polynomial

P (-n,-m) d (x -α) + (x -β) β -α = (2) d j=0 n -d + j -1 j m -j -1 d -j (x -α) j (x -β) d-j (α -β) d ,
and note that it belongs to 1 (α-β) d Z[x -α, x -β] when we consider α and β as distinct indeterminates over Z. We denote by p d its coefficient of x d , for which we show in (14) below that

p d = 1 (α -β) d m + n -d -1 d . (3)
We also recall that, following the notation in Theorem 1, the principal subresultant

s d := PSres d ((x -α) m , (x -β) n ) is the coefficient of x d in Sres d ((x -α) m , (x -β) n ), which by [BDKSV2017, Proposition 3.3] satisfies s d = (α -β) (m-d)(n-d) d i=1 r(i) with r(i) := (i -1)!(m + n -d -i)! (m -i)!(n -i)! . (4)
As a consequence,

s d belongs to Q[α -β] ∩ Z[α, β] = Z[α -β]. In fact it is shown in [BDKSV2017, Theorem 1.1] that the whole polynomial Sres d ((x - α) m , (x -β) n ) belongs to Z[x -α, x -β]
(see also Lemma 6 below for an independent proof). Denote by

Q (-n,-m) d
the following polynomial

Q (-n,-m) d (α, β, x) := s d • P (-n,-m) d (x -α) + (x -β) β -α p d . Since α -β = (x -β) -(x -α), the polynomial Q (-n,-m) d (α, β, x) belongs a priori to Q[x -α, x -β]. We will show that Q (-n,-m) d (α, β, x) actually coincides with Sres d ((x -α) m , (x -β) n ) in Z[x -α, x -β]
and we then obtain, via the map 1 Z → 1 K , the following result: Theorem 2. Let K be a field and α, β ∈ K with α = β. Set d, m, n ∈ N with 0 ≤ d < min{m, n}. Then, with the notation in (2) and (4),

Sres d ((x -α) m , (x -β) n ) = Q (-n,-m) d (α, β, x) .
(5)

The key ingredient to prove Theorem 1 will be to derive from Theorem 2 a second-order recurrence satisfied by the coefficients of Sres d ((x-α) m , (x-β) n ) in the monomial basis, as follows: Theorem 3. Let K be a field and α, β ∈ K with α = β. Set d, m, n ∈ N with 0 ≤ d < min{m, n} and let

Sres d ((x -α) m , (x -β) n ) = d k=0 s k x k .
Then, when char(K) = 0 or char(K) ≥ m + n -d, for s d+1 := 0 and for s d as defined in (4), the following second-order linear recurrence is satisfied by the coefficients s k , for k = d -1, . . . , 0:

s k = -(k + 1) (n -k -1)α + (m -k -1)β s k+1 + (k + 2)αβs k+2 (d -k)(m + n -d -k -1) . (6)
Our next result concerns the complexity of the computation of all principal subresultants PSres d ((x -α) m , (x -β) n ) for 0 ≤ d < min{m, n}. We note that the proof of this result is independent from our previous results, as it is a consequence of a recurrence that is derived directly from (4). We give it here for sake of completeness of our complexity results.

Theorem 4. Let K be a field, let m, n ∈ N and assume char(K) = 0 or char(K) ≥ m + n. Let α, β ∈ K. Then one can compute all the principal subresultants PSres d ((x -α) m , (x -β) n ) ∈ K for 0 ≤ d < min{m, n} using O(min{m, n} + log(mn)) operations in K.

In the current article, we repeatedly use the crucial fact that, for structured algebraic objects, one can obtain improved complexity results by using recurrence relations that these objects obey, rather than just computing them independently. This is one of the strength of our results: not only they provide nice formulae for the subresultants, but they also exploit their particular structure in order to design efficient algorithms.

This work has an interesting story. While working on the paper [BDKSV2017], we first realized that [BDKSV2017, Theorems 1.1 and 1.2] (see Theorem 12 below) implies the linear recurrence on the coefficients of Sres d ((x-α) m , (xβ) n ) in the usual monomial basis described in Theorem 3. This recurrence was initially found using a computer-driven "guess-and-prove" approach, where the guessing part relied on algorithmic Hermite-Padé approximation [START_REF] Salvy | GFUN: a Maple package for the manipulation of generating and holonomic functions in one variable[END_REF], and where the proving part relied on Zeilberger's creative telescoping algorithm [Zei90, WZ1992]. From this we derived a first proof of our complexity result (Theorem 1). Shortly after that, by studying the differential equation attached to this recurrence, we realized that it has a basis of solutions of hypergeometric polynomials, which appeared to be Jacobi polynomials. We have then obtained an indirect and quite involved proof of Theorem 2 and of Theorem 3 based on manipulations of hypergeometric functions, notably on the Chu-Vandermonde identity, much inspired by an experimental mathematics approach. The proof that we choose to present in this article is the shortest and the simplest that we could find. It is chronologically the latest proof of our results, and the one which provides the deepest structural insight. This proof was obtained by applying some classical results and the fact that any polynomial that can be written as a polynomial combination of f and g in K[x] with given degree bounds is in fact a constant multiple of the subresultant of f and g: we prove that the Jacobi polynomial can indeed be expressed as such a combination of (x-α) m and (x -β) n , and we determine the scalar multiple that gives the subresultant. To conclude this introduction, we want to stress here the importance of the interaction between computer science and classical mathematics, which allowed us to guess and prove all our statements using the computer, before finding a short and elegant human proof.

The paper is organized as follows: We first derive Theorems 2 and 3 in Section 2. Section 3 is dedicated to the proof of Theorem 1, while in Section 4 we prove Theorem 4. Section 5 explains the connection of our results with previous work, notably the relationship with classical results on Padé approximation. We conclude the paper with various remarks, experimental results and perspectives in Section 6.

A preliminary version of this work is part of the doctoral thesis of Marcelo Valdettaro [Val2017].

Proofs of Theorem 2 and Theorem 3

Proof of Theorem 2.

The proof of Theorem 2 proceeds in 3 steps: (1) We prove the theorem in the case when K has characteristic 0. (2) We show, independently from [BDKSV2017], that Sres

d ((x -α) m , (x -β) n ) belongs to Z[x -α, x -β] when we consider both polynomials (x -α) m and (x -β) n in Z[α, β, x] for α, β new indeterminates over Z, which implies that (α -β) d p d divides s d • (α -β) d P (-n,-m) d (x -α) + (x -β) β -α in Z[x -α, x -β].
(Here we multiply both terms by (α -β) d to guarantee that they are both polynomials in Z[x -α, x -β].) (3) We finally conclude that the identity stated in Theorem 2 holds in any characteristic via the map 1 Z → 1 K .

We will need the next classical lemma, which follows e.g. from [Mis1993, Lemmas 7.7.4 and 7.7.6] and was also a key ingredient in [BDKSV2017].

Lemma 5. Let m, n ∈ N and f, g ∈ K[x] of degrees m and n respectively. Set 0 ≤ d < min{m, n} and assume Sres d (f, g) = 0 has degree exactly

d. If F, G ∈ K[x] with deg(F) < n-d, deg(G) < m-d are such that h = F f +G g is a non-zero polynomial in K[x] of degree at most d, then there exists λ ∈ K \ {0} satisfying h = λ • Sres d (f, g).
2.1.1. Proof of Theorem 2 when char(K) = 0. In this case Sres d ((x -α) m , (x -β) n ) has degree exactly d by Identity (4) since α = β. We will then show that h = P

(-n,-m) d 2x -α -β β -α satisfies the conditions of Lemma 5 applied to f = (x -α) m and g = (x -β) n .
One can check (or refer to [Sze1975, Theorem 4.23.1] to verify) that the polynomials

P (-n,-m) d (z), (1 + z) m P (-n,m) n-d-1 (z) and (1 -z) n P (n,-m) m-d-1 (z),
all solve the linear differential equation

(1 -z 2 )y (z) + (m + n -2)z -m + n y (z) + d(d + 1 -m -n)y(z) = 0. Substituting z = 2x -α -β β -α in this differential equation shows that the polynomials y 1 (x) := P (-n,-m) d 2x -α -β β -α , y 2 (x) := 2 β -α m (x -α) m P (-n,m) n-d-1 2x -α -β β -α and
y 3 (x) := 2 α -β n (x -β) n P (n,-m) m-d-1 2x -α -β β -α ,
all solve the linear differential equation

(x -α)(x -β)y (x)+ α(n -1) + β(m -1) -(m + n -2)x y (x) + d(m + n -d -1)y(x) = 0. ( 7
)
Since the dimension of the solution space of this second-order linear differential equation is 2, the three polynomials y 1 , y 2 , y 3 must be linearly dependent over K. Now, it is well-known that the Jacobi polynomials satisfy

P (k, ) r (1) = (k + 1) r r! and P (k, ) r (-1) = (-1) r ( + 1) r r! . ( 8
)
This implies that y 2 and y 3 are not linearly dependent over K since

y 2 (β) = 2 m P (-n,m) n-d-1 (1) = (-1) n-d-1 2 m n -1 d = 0 and y 2 (α) = 0, ( 9
)
while

y 3 (β) = 0 and y 3 (α) = 2 n P (n,-m) m-d-1 (-1) = 2 n m -1 d = 0. ( 10
)
Thus, there exist A, B ∈ K such that y 1 (x) = A y 2 (x) + B y 3 (x), that is,

P (-n,-m) d 2x -α -β β -α =A 2 β -α m P (-n,m) n-d-1 2x -α -β β -α (x -α) m (11) + B 2 α -β n P (n,-m) m-d-1 2x -α -β β -α (x -β) n .
In addition P

(-n,-m) d 2x -α -β β -α = 0, since P (-n,-m) d (1) = (-1) d n -1 d and P (-n,-m) d (-1) = m -1 d . (12)
Moreover, deg P

(-n,-m) d 2x-α-β β-α ≤ d, deg P (-n,m) n-d-1 2x-α-β β-α < n -d and deg P (n,-m) m-d-1 2x-α-β β-α < m -d.
Therefore Lemma 5 implies that there exists λ ∈ K such that

P (-n,-m) d 2x -α -β β -α = λ • Sres d ((x -α) m , (x -β) n ). ( 13
)
Thus, the left-hand side and right-hand side of this equality have the same coefficient of x d , which implies that λ = p d /s d . We now determine p d .

By Identity (2),

p d = 1 (α -β) d d j=0 n -d + j -1 j m -j -1 d -j = 1 (α -β) d m + n -d -1 d , (14)
where the second equation can be checked by thinking of a d-combination with repetition from a set of size m + n -2d, written as a disjoint union of a subset with n -d elements and its complement with m -d elements, computed by adding, for 0 ≤ j ≤ d, the j-combination with repetition from the first subset of size n -d combined with the (d -j)-combination with repetition from the second subset of size m -d.

Passing λ -1 = s d /p d to the left-hand side in Identity (13) proves Theorem 2 when char(K) = 0. 2

2.1.2. Proof that Sres d ((x -α) m , (x -β) n ) belongs to Z[x -α, x -β].
This result is already proved in [BDKSV2017], but we give here an independent proof because in Section 5.1 we will show the result in [BDKSV2017] (see Theorem 12 below) and our Theorem 2 are equivalent. Lemma 6. Set d, m, n ∈ N with 0 ≤ d < min{m, n}, and let

(x -α) m , (x - β) n ∈ Z[α, β, x]. Then Sres d ((x -α) m , (x -β) n ) ∈ Z[x -α, x -β].
Proof. It is well-known from the matrix formulation of the subresultant that Sres

d ((x -α) m , (x -β) n ) ∈ Z[α, β, x]. Theorem 2 gives us a way of writing Sres d ((x -α) m , (x -β) n ) = (α -β) (m-d)(n-d) d j=0 c j (x -α) j (x -β) d-j
where c j ∈ Q.

In particular, for α = 0 and β = -1, one has on the one hand

Sres d (x m , (x + 1) n ) = d j=0 c j x j (x + 1) d-j , with c j ∈ Q while on the other hand Sres d (x m , (x + 1) n ) = d k=0 a k x k with a k ∈ Z, 0 ≤ k ≤ d. This means that d j=0 c j x j (x + 1) d-j = d k=0 a k x k , with a k ∈ Z for 0 ≤ k ≤ d.
Comparing coefficients, we deduce that

a k = k j=0 d k -j c j , 0 ≤ k ≤ d, i.e., that    a 0 . . . a d    =      1 d 1 1 . . . . . . . . . d d d d-1 . . . 1         c 0 . . . c d    .
We conclude that c j ∈ Z for all 0 ≤ j ≤ d, since the a k 's are integer numbers and the transition matrix is an invertible integer matrix. 2 2.1.3. Concluding the proof of Theorem 2.

We assume that α and β are distinct indeterminates over Q. The theorem holds over the field Q(α, β), with both sides of equality (5) belonging to Z[x -α, x -β]. To prove the theorem for an arbitrary field K, and for distinct values α and β in K, we apply a classical specialization argument, using the ring homomorphism Z

[x -α, x -β] → K[x] which maps 1 Z → 1 K , α → α, β → β.

Beyond Theorem 2

An advantage of our proof of Theorem 2 is that it also shows that the unique polynomials F d and G d in K[x] of degrees respectively less than n -d and m -d that are the coefficients of the Bézout identity

Sres d ((x -α) m , (x -β) n ) = F d • (x -α) m + G d • (x -β) n , (15)
are also (scalar multiples of) Jacobi polynomials, up to the same affine change of variables. More precisely, we have:

Corollary 7. Let K be a field and α, β ∈ K with α = β. Set d, m, n ∈ N with 0 ≤ d < min{m, n}. Then, the polynomials F d and G d defined in (15) satisfy

F d = (-1) n-1 s d P (-n,m) n-d-1 (x-α)+(x-β) β-α (β -α) m p d , G d = (-1) n s d P (n,-m) m-d-1 (x-α)+(x-β) β-α (β -α) n p d .
Proof. As in the proof of Theorem 2 we first assume that K is a field of characteristic 0. By this theorem, Identities (15) and ( 11), one has

p d F d = s d A 2 β -α m P (-n,m) n-d-1 2x -α -β β -α , p d G d = s d B 2 α -β n P (n,-m) m-d-1 2x -α -β β -α .
We now determine the values of A and B. By Identities (9), ( 10), ( 11) and ( 12), we get

m -1 d = P (-n,-m) d (-1) = B 2 α -β n P (n,-m) m-d-1 (-1)(α -β) n = 2 n m -1 d B, (-1) d n -1 d = P (-n,-m) d (1) = A 2 β -α m P (-n,m) n-d-1 (1)(β -α) m = (-1) n-d-1 2 m n -1 d A.
Therefore A = (-1) n-1 2 m and B = 1 2 n . This proves the statement when char(K) = 0. Finally, both sides in the equalities of Corollary 7 belong to

1 (α-β) m+n-d-1 Z[α, β, x]
and so they specialize well to a field of any characteristic via the map 1 → 1 K . 2 2.3. Proof of Theorem 3. We now prove Theorem 3, which gives a recurrence satisfied by the coefficients (in the monomial basis) of Sres d ((x-α) m , (x-β) n ). The recurrence is inherited from the differential equation (7) satisfied by P

(-n,-m) d (x -α) + (x -β) β -α in characteristic 0.
By Theorem 2,

Sres d ((x -α) m , (x -β) n ) = Q (-n,-m) d (α, β, x) = s d p d • P (-n,-m) d (x -α) + (x -β) β -α , (16)
where

P (-n,-m) d (x -α) + (x -β) β -α
is the integer Jacobi polynomial described in Identity (2), and

s d p d = (α -β) (m-d)(n-d)+d d i=1 r(i) m+n-d-1 d = (α -β) (m-d)(n-d)+d d i=1 i!(m + n -d -i -1)! (m -i)!(n -i)! . (17)
Therefore, the differential equation ( 7) satisfied by the Jacobi polynomial is also satisfied by s(x) := Sres d ((x -α) m , (x -β) n ). We now show that this fact implies the statement. We start with

s(x) = d k=0 s k x k , s (x) = d k=1 ks k x k-1 and s (x) = d k=2 k(k -1)s k x k-2 .
We then have

(x -α)(x -β)s (x) = d k=2 k(k -1)s k x k -(α + β) d k=2 k(k -1)s k x k-1 + αβ d k=2 k(k -1)s k x k-2 = d k=0 k(k -1)s k x k -(α + β) d-1 k=0 (k + 1)ks k+1 x k + αβ d-2 k=0 (k + 2)(k + 1)s k+2 x k , (α(n -1) + β(m -1) -(m + n -2)x) s (x) = -(m + n -2) d k=1 ks k x k + (α(n -1) + β(m -1)) d k=1 ks k x k-1 = -(m + n -2) d k=0 ks k x k + (α(n -1) + β(m -1)) d-1 k=0 (k + 1)s k+1 x k , and
d(m + n -d -1)s(x) = d(m + n -d -1) d k=0 s k x k .
Now we compare the degree-k coefficient in (7) for k = 0, . . . , d -1:

k(k -1) -(m + n -2)k + d(m + n -d -1 s k + -(α + β)(k + 1)k + (α(n -1) + β(m -1))(k + 1) s k+1 + αβ(k + 2)(k + 1)s k+2 = 0.
Therefore,

s k = -(k + 1) (n -k -1)α + (m -k -1)β s k+1 + (k + 2)αβs k+2 (d -k)(m + n -d -k -1) .
This proves the recurrence when char(K) = 0. It is clear that the same recurrence also holds for fields K of characteristic ≥ m + n -d via the map 1 Z → 1 K since in all the steps we are dividing only by natural numbers less than m + n -d. 2

3. Proof of Theorem 1.

Proof of Theorem 1 (a).

We start with the following simple observation.

Lemma 8. Let K be a field, let k, ≥ 0 be integers and assume char(K) = 0 or char(K) > min{k, }. Then the (image in K of the) binomial coefficient 2

The proof that one can compute all coefficients of the d-th subresultant of (x -α) m and (x -β) n in O(min{m, n} + log(mn)) operations in K when char(K) is either zero or larger than m + n -d will be derived from the recurrence (6) described in Theorem 3. The proof is algorithmic and proceeds in several steps. We start with s d = (α -β) (m-d)(n-d) d i=1 r(i), with r(i) defined in (4), and observe that for the mentioned characteristics, s d = 0 since α = β.

• The term (α -β) (m-d)(n-d) can be computed in O(log(mn)) arithmetic operations, by using binary powering.

• The element r(d) = (d-1)! m+n-2d m-d can be computed in O(min{m, n}) arithmetic operations by applying Lemma 8, and using that d < min{m, n}.

• Thanks to the recurrence Note that during the unrolling of the recurrence, the only divisions that occur are by positive integers less than max{m, n}, legitimate in K by the assumption on its characteristic. This shows that s d can be computed using O(min{m, n} + log(mn)) arithmetic operations in K.

r(i) = (m + n -d -i) i(m -i)(n -i) r(i + 1), all r(d -1), . . . , r (
• Starting from s d+1 = 0 and s d , we use the recurrence (6) to compute s d-1 , s d-2 , . . . , s 0 in O(d) operations, by adding O(1) operations in K for each of these d terms.

Note that in this step only divisions by integers less than m + n -d -1 may occur, and all these elements are invertible in K, by assumption.

In conclusion, all the coefficients s 0 , . . . , s d of Sres d ((x -α) m , (x -β) n ) can be computed in O(min{m, n} + log(mn)) operations in K, when char(K) = 0 or char(K) ≥ m + n -d. 2

Proof of Theorem 1 (b).

We apply again the recurrence given by Theorem 3, in the characteristic 0 case, to show that when char(K) = m+n-d-1, the polynomial subresultant Sres d ((x -α) m , (x -β) n )) is actually a (non-zero) constant in K.

Lemma 9. Set d, m, n ∈ N with 1 ≤ d < min{m, n} and let Via the map 1 Z → 1 K , we immediately deduce that s d = • • • = s 1 = 0 in K, and therefore Sres d ((x -α) m , (x -β) n ) ∈ K. We compute its value by specializing Identity (16) at x = α, and thanks to (12) and (17). Set

Sres d ((x -α) m , (x -β) n ) = d k=0 s k x k ∈ Z[α, β][x]. Assume that m + n -d -1 equals a prime number p. Then p | s k in Z[α, β] for 1 ≤ k ≤ d.
p := m + n -d -1 = char(K), then Sres d ((x -α) m , (x -β) n ) is equal to (α -β) (m-d)(n-d)+d d i=1 i!(m + n -d -i -1)! (m -i)!(n -i)! P (-n,-m) d (-1) = (α -β) (m-d)(n-d)+d d i=1 i!(p -i)! (m -i)!(n -i)! m -1 d = (α -β) (m-d)(n-d)+d d i=1 (i -1)!(p -i)! (m -i -1)!(n + i -d -1)! = (α -β) (m-d)(n-d)+d d i=1 p-1 m-i-1 p-1 i-1
.

It remains to show that the last product is equal to (-1) md in K. This is an immediate consequence of the following elementary lemma.

Lemma 10.

p -1 = (-1) in K, for any 0 ≤ < p = char(K).

Proof. By Fermat's little theorem we have (x -1)

p-1 = (x -1) p /(x -1) = (x p -1)/(x -1) = x p-1 + • • • + 1 in K[x]. Thus, the coefficient (-1) p-1 of x in (x -1) p-1 is equal to 1 in K. 2
Finally, by the previous lemma, the following equalities hold in K:

d i=1 p-1 m-i-1 p-1 i-1 = d i=1 (-1) m-i-1 (-1) i-1 = (-1) md .
This concludes the proof of Theorem 1 (b). 2

Proof of Theorem 1 (c).

This non-obvious fact follows for instance from Theorem 2. We know by Identity (16) in the characteristic 0 case that

Sres d ((x -α) m , (x -β) n ) = s d p d • P (-n,-m) d (x -α) + (x -β) β -α where P (-n,-m) d (x -α) + (x -β) β -α
is the integer polynomial described in Identity (2), and

s d p d = (α -β) (m-d)(n-d)+d d i=1 i!(m + n -d -i -1)! (m -i)!(n -i)! .
We note that the denominator in this last term does not vanish in the mentioned characteristics while the numerator equals 0, since it is a multiple of (m + n -d -2)! for d ≥ 1. We conclude the proof of Theorem 1(c) via the map 1 Z → 1 K . 2

Remark. Notice that Theorem 1(c) also follows from Theorem 1(b) and from Collins' fundamental theorem of subresultants ([Col1973, §4], see also [Hab1948, §2] and [Col1967, Theorem 1]) which states that for an arbitrary field K and arbitrary f, g ∈ K[x], the subresultants and the Euclidean remainder sequence of f and g are closely related: if A 1 := f, A 2 := g, A 3 , . . . , A is an Euclidean polynomial remainder sequence of f and g with deg

(A k ) = n k for 1 ≤ k ≤ , then there exist c 1 , . . . , c , d 1 , . . . , d ∈ K × such that Sres n k (f, g) = c k • A k , Sres n k-1 -1 (f, g) = d k • A k , and Sres d (f, g) = 0 for n k < d < n k-1 -1,
for all 1 ≤ k ≤ . In particular, if two nonzero subresultants Sres e (f, g) and Sres e (f, g) have the same degree for some e < e, then they are constant multiples of each other, and all the intermediate subresultants Sres d (f, g) are zero for e < d < e. In our situation, with max{m, n} ≤ p := char(K) < m+n-d-1, and f = (x-α) m , g = (x-β) n in K[x] with α = β, we have that Sres 0 (f, g) ∈ K × and also, by Theorem 1(b), that Sres m+n-p-1 (f, g) ∈ K × . Therefore, Sres d (f, g) = 0 for 1 ≤ d < m + n -1 -p, which reproves part (c) of Theorem 1.

Proof of Theorem 4

With the notation r(i

) := (i -1)!(m + n -d -i)! (m -i)!(n -i)!
introduced in (4), we have:

PSres d ((x -α) m , (x -β) n ) = (α -β) (m-d)(n-d) d i=1 r(i).
While in previous sections d was considered as a fixed value, in this section we view it as variable. Therefore, in order to avoid confusion, we write r d (i) := r(i), to emphasize also its dependence on d. For all integers d ≥ 1, we define

c(d) := d i=1 r d (i)
and note that it is an integer number, as mentioned in the introduction, although the terms r d (i) are not all integers. We also set c(0) := 1.

The key observation for what follows is contained in the next lemma.

Lemma 11. Let K be a field with char(K) = 0 or char(K

) ≥ m + n. Set u(d) := c(d)/c(d -1) for 1 ≤ d < min{m, n} and v(d) := u(d + 1)/u(d) for 1 ≤ d ≤ min{m, n} -2. Then, for 1 ≤ d ≤ min{m, n} -2, v(d) = d(m -d)(n -d)(m + n -d) (m + n -2d -1)(m + n -2d) 2 (m + n -2d + 1) . (18)
At this stage, it remains to compute all the powers h(d) := (α -β) (m-d)(n-d) for 0 ≤ d < min{m, n}, and finally to output PSres

d ((x -α) m , (x -β) n ) = c(d) • h(d)
, for 0 ≤ d < min{m, n}. This is done as follows.

• First, all the elements γ(d) := (α -β) 2d+1-m-n , for d < min{m, n}, are computed using O(log(m + n) + min{m, n}) operations in K. This can be done by computing γ(0) := (α -β) 1-m-n by binary powering, then unrolling the recurrence γ(d+1

) := (α-β) 2 •γ(d) for d < min{m, n}-1.
• Next, h(0) := (α -β) mn is computed by binary powering, and then all h(d), for 1 ≤ d < min{m, n}, by repeated products using h(d + 1) := γ(d) • h(d), for a total cost of O(log(mn) + min{m, n}) operations in K.

• Finally, we compute and return the values PSres

d ((x-α) m , (x -β) n ) = c(d)•h(d), for 0 ≤ d < min{m, n}, using O(min{m, n}) operations in K.
Adding up the various arithmetic costs proves Theorem 4. 2

Connections to previous results

Theorem 2 is closely connected to some previous results. First we discuss the connection to the work [BDKSV2017]. Second, we explain the relationship of the present work to classical results on Padé approximation.

Connection with [BDKSV2017]

We show that the expression for the subresultant obtained in [BDKSV2017], though not expressed in terms of Jacobi polynomials, is equivalent to the one in Theorem 2. First, let us recall the main results of [BDKSV2017].

Theorem 12. [BDKSV2017, Theorems 1.1 and 1.2] Let K be a field and α, β

∈ K. Set d, m, n ∈ N with 0 ≤ d < min{m, n}. Then, Sres d ((x -α) m , (x -β) n ) = (α -β) (m-d)(n-d) d j=0 c j (m, n, d)(x -α) j (x -β) d-j ,
where the coefficients c 0 (m, n, d), . . . , c d (m, n, d) are defined by

c 0 (m, n, d) = d i=1 (i -1)! (m + n -d -i -1)! (m -i -1)!(n -i)! , and
c j (m, n, d) = d j n-d+j-1 j m-1 j c 0 (m, n, d), for 1 ≤ j ≤ d.
(Here c 0 (m, n, 0) = 1, following the convention that an empty product equals 1.) Moreover, for 0

≤ j ≤ d, c j (m, n, d) ∈ Z or Z/pZ if char(K) = 0 or char(K) = p, respectively.
Proof that Theorems 12 and 2 are equivalent. We want to prove that

(α-β) (m-d)(n-d) d j=0 c j (m, n, d)(x-α) j (x-β) d-j = s d P (-n,-m) d 2x -α -β β -α p d , (19)
where

c j (m, n, d) = d j n-d+j-1 j m-1 j d i=1 (i -1)! (c -i)! (m -i -1)!(n -i)! for c := m + n -d -1.
By (17) the right-hand side of ( 19) equals

(α -β) (m-d)(n-d)+d d i=1 i!(c -i)! (m -i)!(n -i)! P (-n,-m) d 2x -α -β β -α ,
where by (2),

(α -β) d P (-n,-m) d 2x -α -β β -α = d j=0 n -d + j -1 j m -j -1 d -j (x -α) j (x -β) d-j .
Thus, we only need to verify that

n -d + j -1 j m -j -1 d -j d i=1 i!(c -i)! (m -i)!(n -i)! = d j n-d+j-1 j m-1 j d i=1 (i -1)!(c -i)! (m -i -1)!(n -i)! ,
i.e. after simplification, that

(m -1)! (m -d-1)! d i=1 i! (m -i)! = d! d i=1 (i -1)! (m -i -1)! ,
which trivially holds. 2

Connection with Padé approximation

In this subsection we show that Theorem 2 and Corollary 7 are also equivalent to classical descriptions of some Padé approximants via Gauss hypergeometric functions.

The starting point is a theorem due to Padé [Pad1901], stating that the

[m/n] Padé approximation in C(x) to (1 -x) k is the ratio of hypergeometric functions 2 F 1 (-m, -k -n; -m -n; x) 2 F 1 (-n, k -m; -m -n; x) . (20)
That result had been previously obtained, by different methods and under several additional assumptions, by Laguerre [Lag1885] and Jacobi [Jac1859].

m (x m+n+1 , (1 -x) k ) G m (x m+n+1 , (1 -x) k ) = (-1) k Sres m (x m+n+1 , (x -1) k ) G m (x m+n+1 , (x -1) k ) , (21)
where

G m := G m (x m+n+1 , (x -1) k ) is the polynomial coefficient of degree ≤ n in the Bézout expression Sres m (x m+n+1 , (x -1) k ) = F m • x m+n+1 + G m • (x -1) k .
Identity (20) implies that

2 F 1 (-m, -n -k; -m -n; x) 2 F 1 (-n, k -m; -m -n; x) = (-1) k Sres m (x m+n+1 , (x -1) k ) G m (x m+n+1 , (x -1) k ) .
We showed earlier that the fact that x m+n+1 and (x-1) k are coprime polynomials implies that deg(Sres m (x m+n+1 , (x-1) k )) = m, and it is also immediate to verify that Sres m (x m+n+1 , (x -1) k ) and G m (x m+n+1 , (x -1) k ) are coprime. Therefore, since the degree of

2 F 1 (-m, -k -n; -m -n; x) = m i=0 (-1) i m i (-k -n) i (-m -n) i x i ,
equals m, one derives that there exists a non-zero λ ∈ C such that

Sres m (x m+n+1 , (x -1) k ) = λ • 2 F 1 (-m, -k -n; -m -n; x), G m (x m+n+1 , (x -1) k ) = (-1) k λ • 2 F 1 (-n, k -m; -m -n; x).
Here, λ can be computed by comparing the leading coefficients of Sres m (x m+n+1 , (x -1) k ) and 2 F 1 (-m, -k -n; -m -n; x):

λ = (-1) m (k + n -m)!(m + n)! (k + n)!n! PSres m (x m+n+1 , (x -1) k ) = (-1) (n+1)(k-m)+m m i=1 (i -1)!(k + n -i)! (k -i)!(m + n -i)! ,
by Identity (4). Now, according to [EMOT1953, (1.6)], see also [Koo1984, (1.5)]:

2 F 1 (-m, -k -n; -m -n; x) = 1 m+n m P (-k,-m-n-1) m (2x -1) , 2 F 1 (-n, k -m; -m -n; x) = 1 m+n m P (k,-m-n-1) n (2x -1) ,
while, according to our Theorem 2 and Corollary 7,

Sres m (x m+n+1 , (x -1) k ) = µ P (-k,-m-n-1) m (2x -1), G m (x m+n+1 , (x -1) k ) = (-1) k µ P (k,-m-n-1) n (2x -1), for µ := (α -β) (m-d)(n-d)+d d i=1 i!(m + n -d -i -1)! (m -i)!(n -i)! and
µ := (-1) (n+1)(k-m)+m m i=1 i!(k + n -i)! (k -i)!(m + n + 1 -i)! .
This shows the equivalence of the results for α = 0, β = 1, since λ = m+n m µ. In order to deduce Theorem 2 and Corollary 7 for any α, β we apply the usual changes of variables formulas that can be found in the now classical book [AJ2006]:

Sres d (f (x -α), g(x -α)) = Sres d (f, g)(x -α), Sres d (f (γx), g(γx)) = γ mn-d(d+1) Sres d (f, g)(γx). Therefore, Sres d ((x -α) m , (x -β) n ) = Sres d (x m , (x -(β -α)) n )(x -α), Sres d (x m , (x -γ) n )(γx) = 1 γ mn-d(d+1) Sres d ((γx) m , (γx -γ) n ) = 1 γ mn-d(d+1) Sres d (γ m x m , γ n (x -1) n ) = γ m(n-d)+n(m-d) γ mn-d(d+1) Sres d (x m , (x -1) n ) = γ (m-d)(n-d)+d Sres d (x m , (x -1) n ).
Hence, since we have just proven that Sres d (x m , (x -

1) n ) = μ P -n,-m d (2x -1) for μ = d i=1 i!(m+n-d-i-1)! (m-i)!(n-i)! , we deduce that Sres d (x m , (x -(β -α)) n )((β -α)x) = μ (β -α) (m-d)(n-d)+d P -n,-m d (2x -1), which implies that Sres d (x m , (x-(β-α)) n )(x) = μ (β-α) (m-d)(n-d)+d P -n,-m d 2 x β -α -1 . We conclude with Sres d ((x -α) m , (x -β) n ) = Sres d (x m , (x -(β -α)) n )(x -α) = μ (β -α) (m-d)(n-d)+d P -n,-m d 2 x -α β -α -1 = μ (β -α) (m-d)(n-d)+d P -n,-m d 2x -α -β β -α ,
as stated in Theorem 2. Note that similar arguments allow to deduce G

d ((x -α) m , (x -β) n ) from G d (x m , (x -1) n ).

Final remarks

Fast computation of cofactors

One can use similar ideas as in the proof of Theorem 1 in order to compute the cofactors F d (x) and G d (x) in Corollary 7 using O(max{m, n} + log(mn)) arithmetic operations in K, when char(K) = 0 or char(K) ≥ max{m, n}. More precisely, we have the following result, whose proof is omitted:

Theorem 13. Let d, m, n ∈ N with 1 ≤ d < min{m, n} and let K be a field with char(K) = 0 or char(K) ≥ max{m, n}, and α, β ∈ K with α = β. Let F d and G d be as defined in (15). Then, (a) if char(K) = 0 or char(K) ≥ m + n -d, then all the coefficients of F d and G d can be computed using O(max{m, n} + log(mn)) arithmetic operations in K, (b) when char(K) = m + n -d -1, the following equalities hold in K

F d = (-1) dm+1 (α -β) (m-d-1)(n-d-1) (x -α) n-d-1 , G d = (-1) dm (α -β) (m-d-1)(n-d-1) (x -β) m-d-1 ,
and the coefficients of F d and G d can be computed using O(max{m, n}+ log(mn)) arithmetic operations in K,

(c) if m + n -d -1 > char(K) ≥ max{m, n} then F d = G d = 0.

Comparison with generic algorithms

As mentioned in the introduction, the fastest algorithms for subresultants of polynomials of degree at most n have arithmetic complexity O(M(n) log n), where M(n) denotes the arithmetic complexity of degree-n polynomial multiplication [Rei1997, LR2001, Lec2018]. These algorithms can compute either one selected polynomial subresultant, or all principal subresultants. Using FFT-based algorithms for polynomial multiplication [GG2013, Ch. 8], their complexity O(M(n) log n) becomes O(n log 2 n log log n), which is quasi-linear up to polylogarithmic factors. These algorithms are generic in the sense that they apply to arbitrary polynomials, and they work in any characteristic.

The algorithms described in the current article are specific to very structured polynomials, namely pure powers of linear polynomials, and they achieve purely linear arithmetic complexity in their maximum degree n. They also compute either one selected polynomial subresultant, or all principal subresultants, but they are restricted to characteristic zero or large enough. The reason is that they require divisions, which is the price to pay for optimality. We leave as an open question whether purely linear arithmetic complexity can be also achieved in arbitrary characteristic.

Another interesting difference is that, while classical algorithms for the order-d subresultant spend more time when d is small (typically, the resultant computation, corresponding to d = 0, is the most expensive), our algorithms spend less time when d is small. For more on practical comparisons, see §6.6.

Algorithmic optimality

The complexity result O(min{m, n} + log(mn)) is quasi-optimal for Theorem 4, since the size of the output is min{m, n}. On the other hand, the complexity result O(min{m, n}+log(mn)) for Theorem 1 is not optimal when d is small compared to m and n. A natural question is whether an algorithm of arithmetic complexity O(d + log(mn)) may exist. While this is true for d = 0, we believe that this is unlikely for d ≥ 1, and moreover we suspect that there is no algorithm for Theorem 1 with arithmetic complexity polynomial in both d and log(mn). Otherwise, we could in particular compute the first principal subresultant PSres 1 ((x -α) m , (x -β) n ) = (α -β) (m-1)(n-1) m + n -2 m -1 , in arithmetic complexity polynomial in log(mn). This does not seem plausible, since it would imply in particular that the central binomial coefficient

2N

N could be computed using an arithmetic complexity polynomial in log N . Although no proof exists, this is generally believed to be impossible.

Fast factorials

It is possible to further improve some of our complexity results by using Strassen's algorithm [Str1976] exploits the special structure of the input polynomials. Table 1 displays some timings for computing Sres d ((x -α) m , (x -β) n ), for various random choices of α, β, m, n and d. Even for moderate degrees m, n, the specialized algorithm is about thousands of times faster. For higher degrees, the generic algorithm becomes quite slow, while the specialized algorithm has a very satisfactory speed.

We also implemented in Maple the algorithm in Section 4, and this time we compared it, on the same examples as in Table 1, with an algorithm written in C by Mohab Safey El Din. The experimental results are displayed in Table 2. Once again, the specialized algorithm is faster than the generic algorithm.

Subresultants for other structured polynomials

The question addressed in this article is a particular case of a much broader topic, the design of efficient algorithms for structured polynomials.

Preliminary results indicate that, for many polynomials whose coefficients satisfy linear recurrences, their subresultants have coefficients that also obey such recurrences; this leaves hope that their computation can be performed in linear time. We plan to study such generalizations in a future work.

  in O(min{k, }) arithmetic operations in K.Proof. It is enough to use for k+ k the most economic of the equivalent writings (k + ) • • • (k + 1)/ ! and ( + k) • • • ( + 1)/k!.

  1) can be deduced from r(d) in O(d) additional operations; then, computing r(1) • • • r(d) also takes O(d) operations.

Proof.

  By applying Identity (4), we first show that p | s d : clearly p does not divide the denominator but p divides (m + n -d -1)! which is in the numerator of r(1). Therefore p | d i=1 r(i) and p | s d (since d ≥ 1). Observe that for 1 ≤ k ≤ d -1, the denominators that appear in the recurrence defining the sequence s k in Theorem 3 range from (d -1)(m + n -d -2) to (m + n -2d), and thus none of them is divisible by p = m + n -d -1. Therefore, since p | s d+1 and p | s d , we inductively conclude that p | s k for 1 ≤ k ≤ d. 2

  See also [Per1913, Eq. (Padé 5), p. 252], [Bak1975, p. 65], [Ise1979] and Theorem 4.1 in [GGZ2012]. There is also a well-known connection between subresultants and Padé approximants (c.f. [GG2013, Corollary 5.21]): the [m/n] Padé approximation in C(x) to (1 -x) k , for integer k ≥ m, equals Sres

  for the computation of N ! in arithmetic complexity O(M( √ N ) log N ), which becomes quasi-linear in √ N when FFTbased algorithms are used for polynomial multiplication. For instance, for fixed d, the principal subresultant PSres d ((x-α) m , (x-β) n ) can be computed using fast factorials in O(d + log(mn) + M( min{m -d, n -d}) log min{m -d, n -d}), operations in K. The same cost can also be achieved for the computation of the whole polynomial subresultant Sres d ((x -α) m , (x -β) n ) in Theorem 1.

Table 2 :

 2 Comparative timings (in seconds) for the computation of the principal subresultants PSres d ((x -α) m , (x -β) n ), on the instances T1-T8 from Table1, using a generic subresultant algorithm implemented in C (column Generic 2), versus the specialized algorithm described in Section 4 implemented in Maple (column New 2). Column Output size G2 displays the bit size of the integer PSres d ((x -α) m , (x -β) n ) computed by Generic 2. Timings displayed in column New 2 correspond to the computation of all PSres k ((x -α) m , (x -β) n ) for 0 ≤ k < min{m, n} -1. Column Output size N2 displays the bit size of the min{m, n} integers computed by New 2.

http://www.regularchains.org/index.html

Acknowledgements. We thank Christian Krattenthaler for precious help with hypergeometric identities during an early stage of this work, and to Mohab Safey El Din for generously sharing his subresultants implementations with us. We are also grateful to the referees for helping us substantially improve the presentation of our results. T. Krick and M. Valdettaro were partially supported by ANPCyT PICT-2013-0294, CONICET PIP-11220130100073CO and UBACyT 2014-2017-20020130100143BA. A. Szanto was partially supported by the NSF grants CCF-1813340 and CCF-1217557.

Proof. We have that u(1) = c(1) = m+n-2 m-1 and for d ≥ 2,

.

Therefore

, which is the desired expression. Note that the only numbers that appear in the denominators of u(d) and of v(d) are products of integers of absolute value less than m + n, which are invertible in K by the assumption on the characteristic of K.

2

Based on Lemma 11, we now design an algorithm that computes all principal subresultants PSres d ((x -α) m , (x -β) n ) with 1 ≤ d < min{m, n} in O(min{m, n} + log(mn)) operations in K, thus proving Theorem 4.

• First, v(1), . . . , v(min{m, n} -2) are computed by using (18) in O(1) arithmetic operations each, for a total of O(min{m, n}) operations in K.

• Then, u(1), . . . , u(min{m, n}-1) are determined, by computing u(1) := m+n-2 m-1

using Lemma 8, in O(min{m, n}) arithmetic operations in K, and by computing iteratively

• Next we compute the elements c(1), . . . , c(min{m, n} -1) iteratively by c(d All examples were run on the same machine, with the latest version of Maple. For entries marked with a -, the computations were aborted after more than 17 hours, with all available memory (150 Gb of RAM) consumed. The last column displays the bit size of the output.

Bit complexity

We have only discussed arithmetic complexity. When K is a finite field, this is perfectly realistic, since arithmetic complexity reflects quite well the running time of the algorithms. When K is infinite, for instance when K = Q, assuming operations in K at unit cost is not realistic anymore, so studying bit complexity becomes a much more pertinent model. Over K = Q, our algorithms in Sections 3 and 4 have very good complexity behaviors in this model too. Indeed, they only involve binary powering, computation of factorials and binomials, unrolling of recurrences, which can be computed in quasi-optimal bit complexity. This is confirmed by the timings in Tables 1 and2, which appear to be indeed quasi-linear in the output size.

Practical issues

The algorithms described in this article have not only a good theoretical complexity, but also a good practical efficiency. We performed some experimental comparisons in Maple, between an implementation of our specialized algorithm in Section 3 and a generic subresultant algorithm available in the package RegularChains 2 . As expected, our algorithm is much faster, since it


J F Bastien

Jens Gustedt

Two's complement sign representation for C2x Modification request for C2x

We implement the agreed change to abandon ones' complement and sign-and-magnitude representation from C

INTRODUCTION

Removing other sign representations than two's complement from C allows to simplify the specification of integer types substantially. As has been voted in the London 2019 meeting we only implement the essential changes to the sign representation in this paper, namely -Remove the specifications of the other sign representations.

-Impose the value of the minimal value of a signed type to be -2 N -1 where N is the width of the type. -Derive the values of the _MIN and _MAX macros for all integer types from the corresponding _WIDTH macro. -Clean-up the remaining parts of the standard from all obsolete mentions of two's complement and negative integer zeros.

WG21 has recently adapted the changes promoted in their document p12363 . Generally, C++ goes much beyond what is presented here:

-Bit-fields can have excess bits.

-Overflowing operations and out-of-range conversion are generally mapped to modulo operations and cannot trap or raise signals.

-Enumeration types and their underlying compatible integer types have precise definitions.

WG14

has not yet found consensus for these points, so we leave them as they are in the current specification.

TWO'S COMPLEMENT

Restricting the possible sign representations to two's complement is relatively straight forward and does not need much of deep thinking. There are some other direct fallouts from doing this, such as other mentions of two's complement in the document that now become obsolete. This concerns in particular the definition of the exact width integer types, and of the (bogus) specifications of arithmetic on atomic types.

TIGHTENING OF INTEGER REPRESENTATIONS

Minimum values of signed integer types

Even for two's complement representation C17 allowed that the value with sign bit 1 and all other bits 0 might be a trap representation. We change this and are thereby in line with the changes in C++. We force that for integer types with a width of N the minimum value is forced to -2 N -1 (and the maximum value remains at 2 N -1 -1).

Adjust widths of signed and unsigned integer types

In C17, the widths of corresponding signed and unsigned may differ by one, in particular an unsigned may be realized by just masking out the sign bit of the signed type. This possibility does not seem to be used in the field, complicates arguing about integers and adds potential case analysis to programs.

WIDTH, MINIMUM AND MAXIMUM MACROS FOR INTEGER TYPES

With these aggreed changes the relationship between the unsigned maximum and signed minimum and maximum values now becomes much simpler and can easily be expressed through the widths of the types, namely if the width is N these values are now fixed to 2 N -1, -2 N -1 and 2 N -1 -1, respectively. Since the integration of the floating point TS's already brings in macros that specify the width of the standard integer types, we simplify the presentation such that it is centered around the width. This has the advantage that all requirements for the minimum width of integer types can now be presented as requirements of _WIDTH macros, and the specification of the _MIN and _MAX can be generic.

PROPOSED TEXT

As usual, we provide a diff-marked set of changed pages in an appendix. Unfortunately, for the central parts of the proposed changes the diff-marked text is not very readable so we provide the whole text for 5.2. 4.2.1, 6.2.6.2, 7.20p5, 7.20.2, and 7.20.3, and only the changes that are textually small are reported in the diffmark appendix.

Characteristics of integer types <limits.h> 1

The values given below shall be replaced by constant expressions suitable for use in #if preprocessing directives. Their implementation-defined values shall be equal or greater to those shown.

-width for an object of type _Bool BOOL_WIDTH 1

-number of bits for smallest object that is not a bit-field (byte)

CHAR_BIT 8

The macros CHAR_WIDTH, SCHAR_WIDTH, and UCHAR_WIDTH that represent the width of the types char, signed char and unsigned char shall expand to the same value as CHAR_BIT. -width for an object of type unsigned short int USHRT_WIDTH 16

The macro SHRT_WIDTH represents the width of the type short int and shall expand to the same value as USHRT_WIDTH. -width for an object of type unsigned int UINT_WIDTH 16

The macro INT_WIDTH represents the width of the type int and shall expand to the same value as UINT_WIDTH. For all unsigned integer types for which <limits.h> or <stdint.h> define a macro with suffix _WIDTH holding its width N , there is a macro with suffix _MAX holding the maximal value 2 N -1 that is representable by the type, that is suitable for use in #if preprocessing directives and that has the same type as would an expression that is an object of the corresponding type converted according to the integer promotions.

3 For all signed integer types for which <limits.h> or <stdint.h> define a macro with suffix _WIDTH holding its width N , there are macros with suffix _MIN and _MAX holding the minimal and maximal values -2 N -1 and 2 N -1 -1 that are representable by the type, that are suitable for use in #if preprocessing directives and that have the same type as would an expression that is an object of the corresponding type converted according to the integer promotions.

4

If an object of type char can hold negative values, the value of CHAR_MIN shall be the same as that of SCHAR_MIN and the value of CHAR_MAX shall be the same as that of SCHAR_MAX. Otherwise, the value of CHAR_MIN shall be 0 and the value of CHAR_MAX shall be the same as that of UCHAR_MAX.4 ...

Integer types 1

For unsigned integer types the bits of the object representation shall be divided into two groups: value bits and padding bits. If there are N value bits, each bit shall represent a different power of 2 between 1 and 2 N -1 , so that objects of that type shall be capable of representing values from 0 to 2 N -1 using a pure binary representation; this shall be known as the value representation. The values of any padding bits are unspecified.The number of value bits N is called the width of the unsigned integer type. There need not be any padding bits; unsigned char shall not have any padding bits.

2

For signed integer types, the bits of the object representation shall be divided into three groups: value bits, padding bits, and the sign bit. If the corresponding unsigned type has width N , the signed type uses the same number of N bits, its width, as value bits and sign bit. N -1 are value bits and the remaining bit is the sign bit. Each bit that is a value bit shall have the same value as the same bit in the object representation of the corresponding unsigned type. If the sign bit N2412:4 JF Bastien and Jens Gustedt is zero, it shall not affect the resulting value. If the sign bit is one, it has value -(2 N -1 ). There need not be any padding bits; signed char shall not have any padding bits.

3

The values of any padding bits are unspecified. A valid (non-trap) object representation of a signed integer type where the sign bit is zero is a valid object representation of the corresponding unsigned type, and shall represent the same value. For any integer type, the object representation where all the bits are zero shall be a representation of the value zero in that type.

4

The precision of an integer type is the number of value bits.

NOTE 1. Some combinations of padding bits might generate trap representations, for example, if one padding bit is a parity bit. Regardless, no arithmetic operation on valid values can generate a trap representation other than as part of an exceptional condition such as an overflow, and this cannot occur with unsigned types. All other combinations of padding bits are alternative object representations of the value specified by the value bits.

NOTE 2. The sign representation defined in this document is called two's complement. Previous revisions of this document additionally allowed other sign representations.

NOTE 3. For unsigned integer types the width and precision are the same, while for signed integer types the width is one greater than the precision.

... For all integer types for which there is a macro with suffix _WIDTH holding the width, maximum macros with suffix _MAX and, for all signed types, minimum macros with suffix _MIN are defined as by 5.2.4.2. ...

Widths of specified-width integer types 1

The following object-like macros specify the width of the types declared in <stdint.h>. Each macro name corresponds to a similar type name in 7.20.1.

2

Each instance of any defined macro shall be replaced by a constant expression suitable for use in #if preprocessing directives. Its implementation-defined value shall be equal to or greater than the value given below, except where stated to be exactly the given value. An implementation shall define only the macros corresponding to those typedef names it actually provides. The following object-like macros specify the width of integer types corresponding to types defined in other standard headers. If it is unspecified if a type is signed or unsigned and the implementation has it as an unsigned type, a minimum macro with extension _MIN, and value 0 of the corresponding type is defined.

2 Each instance of these macros shall be replaced by a constant expression suitable for use in #if preprocessing directives. Its implementation-defined value shall be equal to or greater than the corresponding value given below. An implementation shall define only the macros corresponding to those typedef names it actually provides. 6 7.20.

  Two's complement sign representation for C2x N2412:3 -width for an object of type unsigned long int ULONG_WIDTH 32 The macro LONG_WIDTH represents the width of the type long int and shall expand to the same value as ULONG_WIDTH. -width for an object of type unsigned long long int ULLONG_WIDTH 64 The macro LLONG_WIDTH represents the width of the type long long int and shall expand to the same value as ULLONG_WIDTH. -maximum number of bytes in a multibyte character, for any supported locale MB_LEN_MAX 1 2

  5

  3.1 Width of ptrdiff_t

	N2412:6	JF Bastien and Jens Gustedt
		7.20.3.4 Width of wchar_t
	1	WCHAR_WIDTH	8
		7.20.3.5 Width of wint_t
	1	PTRDIFF_WIDTH	17
		7.20.3.2 Width of sig_atomic_t
	1	SIG_ATOMIC_WIDTH	8
		7.20.3.3 Width of size_t
	1	SIZE_WIDTH	16
	6 A freestanding implementation need not provide all of these types.

1

WINT_WIDTH

16

http://www.open-std.org/jtc1/sc22/wg14/www/docs/n2218.htm

http://www.open-std.org/jtc1/sc22/wg14/www/docs/n2330.pdf

http://wg21.link/p1236 © 2019 by the author(s). Distributed under a Creative Commons Attribution

4.0 International License

See 6.2.5.

Appendix: pages with diffmarks of the proposed changes

The following page numbers are from the particular snapshot and may vary once the changes are integrated.