Factoring sparse polynomials fast

where

, and

for any

. We call

the size of

and

its support. The aim of this paper is to factor

into a product of irreducible sparse polynomials.

All algorithms that we will present are based on the approach of sparse evaluation and interpolation. Instead of directly working with sparse representations (1.1), the idea is to evaluate input polynomials at a sequence of well-chosen points, do the actual work on these evaluations, and then recover the output polynomials using sparse interpolation. The evaluation-interpolation approach leads to very efficient algorithms for many tasks, such as multiplication [47, 45], division, gcd computations [51], etc. In this paper, we investigate the complexity of factorization under this light.

One particularly good way to choose the evaluation points is to take them in a geometric progression: for a fixed

, we evaluate at

, where

. This idea goes back to Prony [88] and was rediscovered, extended, and popularized by Ben Or and Tiwari [5]. We refer to [89] for a nice survey. If

is a finite field, then a further refinement is to use suitable roots of unity, in which case both sparse evaluation and interpolation essentially reduce to discrete Fourier transforms [48, 45].

In this paper, we do not specify the precise algorithms that should be used for sparse evaluation and interpolation, but we will always assume that the evaluation points form geometric progressions. Then the cost

of sparse evaluation or interpolation at

such points is quasi-linear in

. We refer to Sections 2.1, 2.3, and 2.4 for more details on this cost function

and the algebraic complexity model that we assume.

One important consequence of relying on geometric progressions is that this constraints the type of factorization algorithms that will be efficient. For instance, several existing methods start with the application of random shifts

for one or more variables

. Since such shifts do not preserve geometric progressions, this is a technique that we must avoid. On the other hand, monomial transformations like

do preserve geometric progressions and we will see how to make use of this fact.

The main goal of this paper is to develop fast algorithms for factoring sparse polynomials under these constraints. Besides the top-level problem of factorization into irreducibles, we also consider several interesting subtasks, such as gcd computations, Hensel lifting, content-free and square-free factorization, and the extraction of multiple roots. While relying on known techniques, we shall show how to conciliate them with the above constraints.

Our complexity bounds are expressed in terms of the maximal size and total degree of the input and output polynomials. In practical applications, total degrees often remain reasonably small, so we typically allow for a polynomial dependence on the total degree times the required number of evaluation/interpolation points. In this particular asymptotic regime, our complexity bounds are very sharp and they improve on the bounds from the existing literature.

Concerning the top-level problem of decomposing sparse polynomials into irreducible factors, we develop two main approaches: a recursive one on the dimension and a more direct one based on simultaneous lifting with respect to all but one variables. We will present precise complexity bounds and examples of particularly difficult cases.

1.2.Overview of univariate and bivariate factorization methods

The factorization of polynomials is a fundamental problem in computer algebra. Since we are relying on sparse interpolation techniques, it is also natural to focus exclusively on randomized algorithms of Monte Carlo type. For some deterministic algorithms, we refer to [58, 98, 69].

Before considering multivariate polynomials, we need an algorithm for factoring univariate polynomials. Throughout this paper, we assume that we have an algorithm for this task (it can be shown that the mere assumption of

being effective is not sufficient [27, 28]).

In practice, we typically have

, or

for some prime

and

. Most basic is the case when

is a finite field, and various efficient probabilistic methods have been proposed for this case. An early such method is due to Berlekamp [6, 7]. A very efficient algorithm that is also convenient to implement is due to Cantor and Zassenhaus [16]. Asymptotically more efficient methods have been developed since [64, 67] as well as specific improvements for the case when

is composite [53]. See also [31, Chapter 14] and [62].

Rational numbers can either been regarded as a subfield of

. For asymptotically efficient algorithms for the approximate factorizations of univariate polynomials over

, we refer to [95, 84, 80]. When reducing a polynomial in

modulo

for a sufficiently large random prime, factorization over

reduces to factorization over

via Hensel lifting [94, 41, 104]. For more general factorization methods over

, we refer to [26, 20, 79, 38, 3].

Now let

and assume that we have an irreducible factorization

with

for

. (In practice, we require that

are known with sufficient precision.) In order to turn this into a factorization over

, we need a way to recombine the factors, e.g. to find the subsets

for which

. Indeed, if

is irreducible in

and

, then

is never irreducible in

and irreducible with probability

for a large random prime

. The first recombination method that runs in polynomial time is due to Lenstra-Lenstra-Lovasz [74]. Subsequently, many improvements and variants of this LLL-algorithm have been developed [43, 4, 82, 83, 91].

The problem of factoring a bivariate polynomial

over

is similar in many regards to factoring polynomials with rational coefficients. Indeed, for a suitable random prime

, we have seen above that the latter problem reduces to univariate factorization over

, Hensel lifting, and factor combination. In a similar way, after factoring

for a sufficiently random

(possibly in an extension field of

, whenever

is a small finite field), we may use Hensel lifting to obtain a factorization in

, and finally recombine the factors. The precise algorithms for factor recombination are slightly different in this context [29, 68, 4, 71] (see also [90, 92] for earlier related work), although they rely on similar ideas.

1.3.Overview of multivariate factorization methods

Hensel lifting naturally generalizes to polynomials in three or more variables

. Many algorithms for multivariate polynomial factorization rely on it [81, 101, 100, 106, 57, 32, 33, 59, 58, 8, 69, 77, 78, 17], as well as many implementations in computer algebra systems. One important property of higher dimensional Hensel lifting is that factor recombination can generally be avoided with high probability, contrary to what we saw for

-adic and bivariate Hensel lifting. This is due to Hilbert and Bertini's irreducibility theorems [42, 10]:

Effective versions of these theorems can be used to directly reduce the factorization problem in dimension

to a bivariate problem (and even to a univariate problem if

, using similar ideas). We refer to [71] and [72, Chapître 6] for a recent presentation of how to do this and to [97, Section 6.1] for the mathematical background.

In order to analyze the computational complexity of factorization, we first need to specify the precise way we represent our polynomials. When using a dense representation (e.g. storing all monomials of total degree

with their (possibly zero) coefficients), Theorem 1.1 allows us to directly reduce to the bivariate case (if

is very small, then one may need to replace

by a suitable algebraic extension). The first polynomial time reduction of this kind was proposed by Kaltofen [57]. More recent bounds that exploit fast dense polynomial arithmetic can be found in [69].

Another popular representation is the “black box representation”, in which case we are only given an algorithm for the evaluation of our polynomial

. Often this algorithm is actually a straight line program (SLP) [14], which corresponds to the “SLP representation”. In these cases, the relevant complexity measure is the length of the SLP or the maximal number of steps that are needed to compute one black box evaluation. Consequently, affine changes of variables (1.2) only marginally increase the program size. This has been exploited in order to derive polynomial time complexity bounds for factoring in this model [60, 65]; see also [30, 17, 18]. Here we stress that the output factors are also given in black box or SLP representation.

Plausibly the most common representation for multivariate polynomials in computer algebra is the sparse one (1.1). Any sparse polynomial naturally gives rise to an SLP of roughly the same size. Sparse interpolation also provides a way to convert in the opposite direction. However, for an SLP

of length

, it takes

time to recover its sparse representation

, where

. A priori, the back and forth conversion

therefore takes quadratic time

. While it is theoretically possible to factor sparse polynomials using the above black box methods, this is suboptimal from a complexity perspective.

Unfortunately, general affine transformations (1.2) destroy sparsity, so additional ideas are needed for the design of efficient factorization methods based on Hilbert-Bertini-like irreducibility theorems. Dedicated algorithms for the sparse model have been developed in [32, 8, 76, 78]. There are two ways to counter the loss of sparsity under affine transformations, both of which will be considered in the present paper. One possibility is to successively use Hensel lifting with respect to

. Another approach is to use a more particular type of changes of variables, like

. However, both approaches require

to be of a suitable form to allow for Hensel lifting. Some reference for Hensel lifting in the context of sparse polynomials are [59, 8, 75, 77, 78, 17].

For most applications in computer algebra, the total degree of large sparse polynomials is typically much smaller than the number of terms. The works in the above references mainly target this asymptotic regime. The factorization of “supersparse” polynomials has also been considered in [39, 36]. The survey talk [93] discusses the complexity of polynomial factorizations for yet other polynomial representations.

The theory of polynomial factorization involves several other basic algorithmic tasks that are interesting for their own sake. We already mentioned the importance of Hensel lifting. Other fundamental operations are gcd computations, multiple root computations, square-free factorizations, and determining the content of a polynomial. We refer to [31] for classical univariate algorithms. As to sparse multivariate polynomials, there exist many approaches for gcd computations [21, 105, 65, 63, 56, 22, 54, 73, 51, 55].

Whenever convenient, we will assume that the characteristic of

is either zero or sufficiently large. This will allow us to avoid technical non-separability issues; we refer to [34, 99, 70] for algorithms to deal with such issues. A survey of multivariate polynomial factorization over finite fields (including small characteristic) can be found in [62]. Throughout this paper, factorizations will be done over the ground field

. Some algorithms for “absolute factorization” over the algebraic closure

can be found in [61, 19, 23]; alternatively, one may mimic computations in

using dynamic evaluation [24, 49].

1.4.Outline of our contributions

The goal of this paper is to redevelop the theory of sparse polynomial factorization, by taking advantage as fully as possible of evaluation-interpolation techniques at geometric sequences. The basic background material is recalled in Section 2. We recall that almost all algorithms in the paper are randomized, of Monte Carlo type. We also note that the correctness of a factorization can easily be verified, either directly or with high probability by evaluating the polynomial and the product of the potential factors at a random point.

As an appetizer, we start in Section 3 with the problem of multivariate gcd computations. This provides a nice introduction for the main two approaches used later in the paper: induction on dimension and direct reduction to the univariate (or bivariate or low dimensional) case. It also illustrates the kind of complexity bounds that we are aiming for. Consider the computation of the gcd

of two sparse polynomials

. Ideally speaking, setting

, we are aiming for a complexity bound of the form

where

is a constant. Since

is typically much smaller than

, we can afford ourselves the non-trivial overhead with respect to

in the term

. The inductive approach on the dimension

achieves the complexity (1.3) with

on many practical instances, but its worst case complexity is

. This algorithm seems to be new. Using a direct reduction to univariate gcd computations via so-called “regularizing weights”, the second approach achieves the complexity (1.3) with

, and even

for many practical examples. This algorithm is a sharpened version of the algorithm from [51, Section 4.3].

Most existing algorithms for multivariate polynomial factorization reduce to the univariate or bivariate case. Direct reduction to the univariate case is only possible for certain types of coefficients, such as integers, rational numbers, or algebraic numbers. Reduction to the bivariate case works in general, thanks to Theorem 1.1, and this is even interesting when

: first reduce modulo a suitable prime

, then factor over

, and finally Hensel lift to obtain a factorization over

. In this paper, we will systematically opt for the bivariate reduction strategy. For this reason, we have included a separate Section 4 on bivariate factorization and related problems. This material is classical, but recalled for self-containedness and convenience of the reader.

is a finite field, then we already noted that multivariate factorization does not directly reduce to the univariate case. Nevertheless, such direct reductions are possible for some special cases of interest: content-free factorization, extraction of multiple roots, and square-free factorization. In Section 5, we present efficient algorithms for these tasks, following the “regularizing weights” approach that was introduced in Section 3.2 for gcd computations. All complexity bounds are of the form (1.3) for the relevant notions of output size

, input-output size

, and total degree

In Section 6 we turn to the factorization of a multivariate polynomial

using induction on the dimension

. Starting with a coprime factorization of a bivariate projection

for random

, we successively Hensel lift this factorization with respect to

. Using a suitable evaluation-interpolation strategy, the actual Hensel lifting can be done on bivariate polynomials. This leads to complexity bounds of the form

with

. In fact, we separately consider factorizations into two or more factors. In the case of two factors, it is possible to first determine the smallest factor and then perform an exact division to obtain the other one. In the complexity bound, this means that

should really be understood as the number of evaluation-interpolation points, i.e. the minimum of the sizes of the two factors. It depends on the situation whether it is faster to lift a full coprime factorization or one factor at a time, although we expect the first approach to be fastest in most cases.

Due to the fact that projections such as

are only very special types of affine transformations, Theorem 1.1 does not apply. Therefore, the direct inductive approach from Section 6 does not systematically lead to a full irreducible factorization of

. In Remarks 6.14 and 6.15, we give an example on which our approach fails, together with two different heuristic remedies (which both lead to similar complexity bounds, but with

or higher).

In the last Section 7, we present another approach, which avoids induction on the dimension

and the corresponding overhead in the complexity bound. The idea is to exploit properties of the Newton polytope of

and “face factorizations” (e.g. factorizations of restrictions of

to faces of its Newton polytope). In the most favorable case, there exists a coprime edge factorization, which can be Hensel lifted into a full factorization, and we obtain a complexity bound of the form (1.3). In less favorable cases, different face factorizations need to be combined. Although this yields a similar complexity bound, the details are more technical. We refer to [1, 102] for a few existing ways to exploit Newton polytopes for polynomial factorization.

In very unfavorable cases (see Section 7.6), the factorization of

cannot be recovered from its face factorizations at all. In Section 7.7, we conclude with a fully general algorithm for irreducible factorization. This algorithm is not as practical, but addresses pathologically difficult cases through the introduction of a few extra variables. Its cost is

plus the cost of one univariate factorization of degree

1.5.Notation

In this paper paper, we will measure the complexity of algorithms using the algebraic complexity model [14]. In addition, we assume that it is possible to sample a random element from

(or a subset of

) in constant time. We will use

as a shorthand for

We let

and

. We also define

, for any ring

. Given a multivariate polynomial

and

, we write

(resp.

) for the degree (resp. valuation) of

. We also write

(resp.

) for the total degree (resp. valuation of

), and we set

. Recall that

stands for the number of terms of

in its sparse representation (1.1).

Acknowledgment. We wish to thank Grégoire Lecerf for useful discussions during the preparation of this paper.

2.Preliminaries

2.1.Basic complexities

Let

(or

) be the cost to multiply two dense univariate polynomials of degree

. Throughout the paper, we make the assumption that

is a non-decreasing function. In the algebraic complexity model [14], when counting the number of operations in

, we may take

[15]. If

is a finite field

, then one has

, under suitable number theoretic assumption [40]. In this case, the corresponding bit complexity (when counting the number of operations on a Turing machine [85]) is

For polynomials

of degree

it is possible to find the unique

, such that

with

. This is the problem of univariate division with remainder, which can be solved in time

by applying Newton iteration [31, Section 9.1]. A related task is univariate root extraction: given

and

, check whether

is of the form

for some

and monic

, and determine

and

if so. For a fixed

, this can be checked, and the unique

can be found in

arithmetic operations in

by applying Newton iteration [31, polynomial analogue of Theorem 9.28].

where

are pairwise distinct. Given a column vector

with entries

, it is well known [11, 66, 12, 9, 44] that the products

, and

can all be computed in time

. These are respectively the problems of multi-point evaluation, (univariate) interpolation, transposed multi-point evaluation, and transposed interpolation. For fixed

, these complexities can often be further reduced to

using techniques from [44].

For our factorization algorithms, we will sometimes need to assume the existence of an algorithm to factor univariate polynomials in

into irreducible factors. We will denote by

the cost of such a factorization as a function of

. We will always assume that

is non-decreasing. In particular

for all

and

is the finite field

with

elements for some odd

, then the best univariate factorization methods are randomized of Las Vegas type. When allowing for such algorithms, we may take

, by using Cantor and Zassenhaus' method from [16]. With the use of fast modular composition [67], we may take

but this algorithm is only relevant in theory, for the moment [50]. If the extension degree

is composite, then this can be exploited to lower the practical complexity of factorization [53].

2.2.The Schwarz–Zippel lemma

In all probabilistic algorithms in this paper,

will stand for a sufficiently large integer such that “random elements in

” are chosen in some fixed subset of

of size at least

. In the case when

is larger than the cardinality

, this means that

needs to be replaced by an algebraic extension of degree

, which induces a logarithmic

overhead for the cost of field operations in

. We will frequently use the following well-known lemma:

2.3.Sparse polynomial interpolation

Consider a polynomial

that is presented as a blackbox function. The task of sparse interpolation is to recover the sparse representation (1.1) from a sufficient number of blackbox evaluations of

. One may distinguish between the cases when the exponents

are already known or not. (Here “known exponents” may be taken liberally to mean a superset of reasonable size of the set of actual exponents.)

One popular approach for sparse interpolation is based on Prony's geometric sequence technique [88, 5]. This approach requires an admissible ratio

, such that for any

, there is an algorithm to recover

from

. If

, then we may take

to be the

-th prime number, and use prime factorization in order to recover

. If

is a finite field, where

is a smooth prime (i.e.

has many small prime divisors), then one may recover exponents using the tangent Graeffe method [37].

Given an admissible ratio

, Prony's method allows for the sparse interpolation of

using

evaluations of

for

, as well as

operations in

for determining

, and

subsequent exponent recoveries. If

, then the exponents can be recovered if

, which can be ensured by working over a field extension

with

. If the exponents

are already known, then the coefficients can be obtained from

evaluations of

using one transposed univariate interpolation of cost

is a finite field, then it can be more efficient to consider an FFT ratio

for which

are roots of unity. When choosing these roots of unity with care, possibly in an extension field of

, sparse interpolation can often be done in time

from

values of

, using discrete Fourier transforms; see [48, 45] for details.

In what follows, we will denote by

the cost of sparse interpolation of size

, given

values of

at a suitable geometric sequence. When using Prony's approach, we may thus take

, whenever the cost to recover the exponents

from

is negligible. If the discrete Fourier approach is applicable, then we may even take

Remark 2.4. If we have a bound

for the number of terms of

, then we assume that our sparse interpolation method is deterministic and that it interpolates

in time

. If we do not know the number of terms of

, then we may run the interpolation method for a guessed number of

terms. We may check the correctness of our guessed interpolation

by verifying that the evaluations of

and

coincide at a randomly chosen point. By the Schwartz-Zippel lemma, this strategy succeeds with probability at least

Remark 2.5. The above interpolation methods readily generalize to the case when we use a geometric progression of the form

with

as the evaluation points, by considering the function

instead of

. Taking a random

avoids certain problems due to degeneracies; this is sometimes called “diversification” [35]. If

is itself chosen at random, then it often suffices to simply take

. We will occasionally do so without further mention.

Remark 2.6. For many probabilistic proofs in the sequel of this paper, we will rely on Corollary 2.2. However, this requires the ratios

of geometric progressions to be picked at random, which excludes FFT ratios. Alternatively, we could have relied on Corollary 2.3 and diversification of all input and output polynomials for a fixed random scaling factor

(see the above remark).

Assume for instance that we wish to factor

. Then the factors

are in one-to-one correspondence with the factors

. If we rely on Corollary 2.3 instead of 2.2 in our complexity bounds for factoring

(like the bounds in Theorems 6.11 or 7.2 below), then the cost of diversification adds an extra term

, where

and

is the total size of

and the computed factors. On the positive side, in the corresponding bounds for the probability of failure, the quadratic dependence

on the number

of evaluation points reduces to a linear one. Similar adjustments apply for other operations such as gcd computations.

2.4.Sparse evaluation at geometric progressions

Opposite to sparse interpolation, one may also consider the evaluation of

points

in a geometric progression. In general, this can be done in time

, using one transposed multi-point evaluation of size

. If

is a suitable FFT ratio, then this complexity drops to

, using a discrete Fourier transform of size

. In what follows, we will assume that this operation can always be done in time

More generally, we may consider the evaluation of

points

in a geometric progression. If

, we may do this in time

, by reducing to the evaluation of

progressions of size

. To obtain the evaluations of

for

, we can evaluate

. If

, then we may cut

into

polynomials of size

, and do the evaluation in time

Remark 2.7. If

is an FFT-ratio, then the bound for the case when

further reduces to

plus

bit operations on exponents, since the cyclic projections from [45, 48] reduce

in linear time to cyclic polynomials with

terms before applying the FFTs. We did not exploit this optimization for the complexity bounds in this paper, since we preferred to state these bounds for general ratios

, but it should be straightforward to adapt them to the specific FFT case.

In this paper, we will frequently need to evaluate

at all but one or two variables. Assume for instance that

where

has

terms for

and

. Then using the above method, we can compute

for

using

operations.

One traditional application of the combination of sparse evaluation and interpolation are probabilistic algorithms for multiplication and exact division of sparse multivariate polynomials. For the given

, we can compute the product

by evaluating

, and

at a geometric progression

with

and recovering of

in the sparse representation (1.1) via sparse interpolation. The total cost of this method is bounded by

operations in

, where

. If

and

are known, then the exact quotient

can be computed in a similar fashion and with the same complexity. If

, then the quotient

can actually be computed in time

. Divisions by zero are avoided through diversification, with probability at least

, by Corollary 2.2.

2.5.Newton polytopes

Consider a multivariate polynomial

. We define

to be the convex hull of

and we call it the Newton polytope of

. Given another polynomial

, it is well known that

Let

be a non-zero weight vector. We define the -degree, -valuation, and -ecart of a non-zero polynomial

Note that

and

, where we exploit the fact that we allow for negative weights. We say that

is -homogeneous if

. Any

can uniquely be written as a sum

We call

and

the -leading and -trailing parts of

. We say that

is -regular if

consists of a single term

. In that case, we denote

and we say that

is -monic if

. Given another non-zero polynomial

, we have

The Newton polytopes

and

are facets of

. In particular, they are contained in the boundary

2.6.Laurent polynomials

Consider the rings

and

of ordinary polynomials and Laurent polynomials. Both rings are unique factorization domains, but the group of units of

, whereas the group of units of

. Factoring in

is therefore essentially the same as factoring in

up to multiplications with monomials in

. For instance, the factorization

gives rise to the factorization

. Conversely, the factorization

gives rise to the factorization

. Similarly, computing gcds in

is essentially the same problem as computing gcds in

This is clearly a homomorphism, which is injective (or surjective) if the linear map

is injective (or surjective). Note also that

for any matrices

and

. In particular, if

is unimodular, then

is an automorphism of

with

Laurent polynomials are only slightly more general than ordinary polynomials and we already noted above that factoring in

is essentially the same problem as factoring in

(and similarly for gcd computations). It is also straightforward to adapt the definitions from Section 2.5 and most algorithms for sparse interpolation to this slightly more general setting. The main advantage of Laurent polynomials is that they are closed under monomial maps, which allows us to change the geometry of the support of a polynomial without changing its properties with respect to factorization.

2.7.Tagging

This map is an injective monomial map. For any

, we have

, and

. Divisibility and gcds are preserved as follows:

Proof. We claim that

divides

if and only if

divides

. One direction is clear. Assume that

divides

, so that

with

. We may uniquely write

with

and

such that

for

. Since

for

, it follows that

. Hence

divides

. Our claim implies that

divides

if and only if

divides

. Now we may further extend

to a monomial automorphism of

by sending

to itself. This yields (a). The second property is an easy consequence.

Given a

-regular polynomial

, we note that any divisor

must again be

-regular. Hence, modulo multiplication with a suitable constant in

, we may always normalize such a divisor to become

-monic. In the setting of Laurent polynomials, we may further multiply

by a monomial in

such that

. Similarly, when applying

, we can always normalize

to be monic as a Laurent polynomial in

by considering

Both for gcd computations and factorization, this raises the question of how to find weights

for which a given non-zero polynomial

-regular. In [51, Section 4.3], a way was described to compute such a regularizing weight

: let

be such that

is maximal. In that case, it suffices to take

, and we note that

, where

is the total degree of

. For our applications in Sections 3.2 and 5 it is often important to find a

for which

is small. By trying a few random weights with small (possibly negative) entries, it is often possible to find a regularizing weight

with

Example 2.10. Consider

. Then,

is not

-regular for the natural weight

. If we take

instead, then

-regular with

, and we arrive at

Furthermore,

can be normalized to be monic as a Laurent polynomial in

by considering

. Note that

is also a regularizing weight for

with

3.Multivariate gcd computations

Before studying the factorization problem for multivariate polynomials, it is interesting to consider the easier problem of gcd computations. In this section we introduce two approaches for gcd computations that will also be useful later for factoring polynomials.

The first approach is iterative on the number of variables. It will be adapted to the factorization problem in Section 6. The second approach is more direct, but requires a regularizing weight (see Section 2.7). Square-free factorization can be accomplished using a similar technique, as we shall see in Section 5. The algorithms from Section 7 also draw some of their inspiration from this technique, but also rely on Newton polygons instead of regularizing weights.

3.1.Iterative computation of gcds

Let

and

. As we will see below,

for

with high probability. We may easily compute the univariate greatest common divisor

. In this subsection, we shall describe an iterative algorithm to compute

from

for

. Eventually, this yields

For any

, we have

with high probability. Now these greatest common divisors are only defined up to non-zero scalar multiples in

. Nevertheless, there exists a unique greatest common divisor

and

whose evaluation at

coincides with

is known, then

, and

can be computed for successive

using fast evaluation at geometric progressions. For any

, we may then compute the univariate gcd

and

, under the normalization constraint that

. It finally suffices to interpolate

from sufficiently many

. This yields the following algorithm:

If , then compute using a univariate algorithm and return it

Compute by recursively applying the algorithm to and

Compute , , for using sparse evaluation

Recover from using sparse interpolation

Before we analyze this algorithm, we will need a few probabilistic lemmas. Assume that

and

are independently chosen at random from a subset of

of size at least

(if

is a small finite field, this forces us to move to a field extension).

Proof. Let us write

and

, where

and

are coprime. Then we have

if and only if

and

are not coprime. The result now follows by applying the previous lemma for

Remark. The probability bound implicitly assumes that

, since the statement becomes void for smaller

. In particular, we recall that this means that the cardinality of

should be at least

Proof. Assuming that

and

for

and

, let us prove that Algorithm 3.1 returns the correct answer. Indeed, these assumptions imply that

can be normalized such that

and that the unique such

must coincide with

. Now

fails for some

with probability at most

, by Lemma 3.2. Since

, the condition

fails with probability at most

for some

and

, by Corollary 2.2. This completes the probabilistic correctness proof.

As to the complexity, let us first ignore the cost of the recursive call. Then the sparse evaluations of

, and

can be done in time

: see Section 2.4. The univariate gcd computations take

operations. Finally, the recovery of

using sparse interpolation can be done in time

. Altogether, the complexity without the recursive calls is bounded by

. We conclude by observing that the degrees and numbers of terms of

, and

can only decrease during recursive calls. Since the recursive depth is

, the complexity bound follows.

Remark 3.4. When recovering

from

using sparse interpolation, one may exploit the fact that the exponents of

are already known.

Example 3.5. Let

be random polynomials of total degree

and consider

. With high probability,

. Let us measure the overhead of recursive calls in Algorithm 3.1 with respect to the number of variables

. With high probability, we have

This shows that the sizes of the supports of the input and output polynomials in Algorithm 3.1 become at least twice as small at every recursive call. Consequently, the overall cost is at most twice the cost of the top-level call, roughly speaking.

3.2.Gcd computations through regularizing weights

Algorithm 3.1 has the disadvantage that the complexity bound in Theorem 3.3 involves a factor

. Let us now present an alternative algorithm that avoids this pitfall, but which may require a non-trivial monomial change of variables. Our method is a variant of the algorithm from [51, Section 4.3].

Given

, we first compute a regularizing weight

for

, for which

is as small as possible. From a practical point of view, as explained in Section 2.7, we first try a few random small weights

. If no regularizing weight is found in this way, then we may always revert to the following choice:

Proof. Assume that

. Then

for any

with

. The case when

is handled similarly.

Now consider

, and

. We normalize

in such a way that

and such that

is monic as a polynomial in

; this is always possible since

-regular. Let

be an admissible ratio or an FFT ratio. For any

, we evaluate at

and

, which leads us to define the univariate polynomials

With high probability,

, and

is the monic gcd of

and

, where

and

. For sufficiently large

(with

), we may thus recover

from

using sparse interpolation. Finally, for

, we obtain

, where

for

. This leads to the following algorithm:

Output: , such that

Find a regularizing weight for or

Compute with monic and for

If yield through sparse interpolation, then

Remark 3.7. If

is normalized in

to be monic as a polynomial in

, then we may need to interpolate multivariate polynomials with negative exponents in order to recover

from

. In practice, many interpolation algorithms based on geometric sequences, like Prony's method, can be adapted to do this.

As in the previous subsection, assume that

are independently chosen at random from a subset of

of size at least

. We let

, and

Proof. We have

if and only if

does not vanish at

. Now the degree of

is at most

, so the probability that

for a randomly chosen

and all

is at least

, by Corollary 2.2.

Proof. The correctness with the announced probability follows from Corollaries 2.9 and 3.9, while also using Remark 2.4. The computation of

and

through sparse evaluation at geometric progressions requires

operations (see Section 2.4). The univariate gcd computations take

further operations. The final interpolation of

from

can be done in time

Example 3.11. Let

. Consider the particular case when

is monic as a polynomial in

. Then,

is a regularizing weight for

, and therefore also for

. This situation can be readily detected and, in this case, we have

in the complexity bound (3.1).

Remark 3.12. We may need fewer evaluation points to interpolate the gcd

in Algorithm 3.2 in case the terms of

are distributed over the powers of

. For instance, if the terms are distributed evenly, then we have

in the complexity bound (3.1).

4.Bivariate factorization

Lemma 3.2 allows us to project the general problem of multivariate gcd computations down to the univariate case. For the polynomial factorization problem, no such reduction exists: given a random univariate polynomial of degree

over a finite field, there is a non-zero probability that this polynomial is not irreducible. For this reason, it is customary to project down to bivariate instead of the univariate polynomials (when applicable, an alternative is to project down to univariate polynomials with integer coefficients; see [76], for instance).

This explains why it is interesting to study the factorization of bivariate polynomials in more detail. Throughout this section,

is a bivariate polynomial in

of degree

and of degree

. As in Sections 2.2 and 3, random numbers will be drawn from a subset of

of size at least

. We will recall some classical results concerning the complexity of bivariate factorization, important techniques, and special cases: content-free factorization, root extraction, Hensel lifting, and square-free factorization.

4.1.Content-free factorization

We say that

is content-free (or primitive) in

. Given two random shifts

, we have

with high probability. More precisely:

Proof. Without loss of generality, we may assume that

. Let us first consider the case when

. Then we claim that

. Indeed, for

pairwise distinct

, the Vandermonde matrix

is invertible, so

. It follows that

, regarded as an element of

, is non-zero, and its total degree is bounded by

. By Lemma 2.1, it follows that

with probability at least

. In that case,

We have proved our probabilistic correctness claim in the particular case when

. In general, we factor

. With probability at least

, we have

As to the complexity bound, the evaluations

and

require

operations and the univariate gcd computation can be done in time

4.2.Root extraction

Let

and

. Assume that

for some

and

. Assume that

. Then

can be computed efficiently as follows.

After dividing out a suitable power of

, we may assume without loss of generality that

. For a random shift

, we next replace

with

. With high probability, this ensures that

. Modulo division of

, we may then assume without loss of generality that

Let

be an admissible ratio or an FFT ratio. For any

, we define the univariate polynomial

. With high probability, we have

. Let

be the unique univariate polynomial such that

. Such

can be found efficiently using univariate root extraction. For

, we may recover

from

using interpolation.

Proof. The random shift

and the corresponding backshift

can be computed in time

using the so-called convolution method [2, Theorem 5]. The

can also be computed in time

using fast multipoint evaluation at geometric sequences [2, 13]. The same holds for the interpolation of

from the

. Finally, the univariate root extractions can be computed in time

. The algorithm is correct as long as

and

for

. The probability that

and

for a random choice of

is at least

, by Lemma 2.1. In that case,

and the probability that

for

is at least

, by Corollary 2.2.

4.3.Hensel lifting

Let

be content-free in

and assume that

has a non-trivial factorization

, and

. Assume that

and that

and

are known and coprime (in particular

and

are coprime). Using a random shift

, these assumptions can be enforced with high probability, provided that

and

are coprime. Without loss of generality we may also assume that we normalized the factorization of

such that

is monic. Under these assumptions, we may compute

and

as follows:

Proof. We first observe that there is a unique factorization

that lifts the factorization of

, thanks to our hypothesis that

. Since this hypothesis also implies that

, the denominator

as an element of

coincides with the leading coefficient of

as a polynomial in

. Consequently, the denominator of

equals

if and only if

. Since the degree of

is bounded by

, this happens with probability at least

. Since the degrees of the numerator

and denominator

do not exceed

, it suffices to compute

modulo

in order to recover

and

. This completes the probabilistic correctness proof.

As to the complexity bound, the Hensel lifting requires

operations in

, when using a fast Newton iteration [31, Theorem 15.18, two factors]. The computation of

requires

further operations and the rational function reconstruction can be done in time

using the technique of half gcds [31, Chapter 11].

The proposition generalizes in a straightforward way to the case when

has

pairwise coprime factors

. In that case, one may use fast multifactor Hensel lifting [31, Theorem 15.18], to obtain the following:

4.4.Square-free factorization

Assume that

is content-free in

and of total degree

. Assume also that

. Recall that the square-free factorization of

is of the form

where each

is the product of all irreducible factors of

that occur with multiplicity

. Note that some of the

are allowed to be units in

and that the

are unique up to multiplication by such units. The polynomials

are pairwise coprime. Since

, they must also be separable in both

and

(i.e.

). The square-free factorization of

can be computed efficiently as follows:

Proof. Given

, consider

and

. The polynomials

and

are coprime if and only if

does not vanish at

. Since this resultant has degree at most

, this happens with probability

. Therefore, the probability that all

are pairwise coprime and all

are separable is at least

. In that case,

is the square-free decomposition of

. Modulo normalization, we are thus in the position to apply Proposition 4.4. This proves the probabilistic correctness of the algorithm.

The computation of the shift

can be done in time

using the algorithm from [2, Theorem 5] and the same holds for the shifts in the opposite direction in the last step. The square-free factorization of the univariate polynomial

can be done in time

: see [103] and [31, Theorem 14.23]. We conclude with Proposition 4.4.

4.5.General bivariate factorization

General bivariate factorization of

can often be reduced to Hensel lifting as in Section 4.3 using a random shift

and diversification

. Let

. The best currently known complexity bound is the following:

The actual proposition in [71] also contains a similar result for finite fields of small characteristic. For

as in Theorem 4.7, square-freeness and content-freeness can be achieved with high probability and negligible cost using the algorithms from Sections 4.1 and 4.4.

5.Efficient reductions

In the bivariate setting of the previous section, we have presented several efficient algorithms for the computation of partial factorizations. In this section, we will generalize three of them to the multivariate setting: removal of content, root extraction, and square-free factorizations. The common feature of these generalizations is that they recover the answer directly from the corresponding univariate specializations of the problem, in a similar fashion as the gcd algorithm from Section 3.2.

5.1.Content-free factorization

Consider a polynomial

and a variable

. If, for every non-trivial factorization

with

, both

and

depend on

, then we say that

is content-free (or primitive) in

. Note that this is always the case if

. If

is content-free in

for all

, then we say that

is content-free.

For a given variable

, say

, we can efficiently test whether

is content-free with respect to

: for random

, we form the bivariate polynomial

and compute

using the method from Section 4.1. With high probability,

is content-free with respect to

if and only if

. Performing this test for each of the variables

yields:

Proof. The proof is similar to the one of Proposition 4.1. This time, for the probability bound, we consider the resultant

as a polynomial in

, of total degree at most

. If

, then this resultant does not vanish with probability at least

for random

As to the complexity bound, for

, let

and

be random and consider

and

. We compute the

simultaneously for

in time

and similarly for the

. Finally, the computation of

for

takes

operations.

Assume now that

is not content-free, say with respect to

. With high probability, the content of

with respect to

then equals the gcd of

and

, for two random shifts

. This leads to a non-trivial factorization of

for the cost of one gcd computation and one exact division.

5.2.Root extraction

Given

and

, multivariate

-th root extraction is the problem of determining

and

, such that

, whenever such

and

exist. We devise an algorithm in the same vein as the algorithm for gcd computations from Section 3.2.

We first compute a regularizing weight

for

such that

is small. Recall that the regularity of

ensures that

for some

and

. We take

and note that we then must have

Now let

with

, so that

and

is monic as a polynomial in

. Let

be an admissible ratio or an FFT ratio. For any

, we define the univariate polynomial

. Let

be the unique monic polynomial with

. For sufficiently large

(with

), we may recover

from

using sparse interpolation. Finally, we have

, where

is such that

for

Proof. The evaluations

take

operations, whereas the sparse interpolation of

from the

can be done in time

. The cost of the univariate

-th root extractions

is bounded by

5.3.Square-free factorization

Consider

of total degree

. Assume that

is content-free and that

. The factorization of

into square-free parts can be done using a similar technique as for gcd computations in Section 3.2. We start with the computation of a regularizing weight

for

. Setting

, we recall that

, whence

. Let

where

, and where

are monic and of valuation zero in

. Let

be an admissible ratio or an FFT ratio. For any

, we define the univariate polynomials

where

, and

are monic polynomials in

. We recover

using sparse interpolation and then

by setting

and multiplying by appropriate monomials in

. More precisely,

, where

and

otherwise.

Proof. The probabilistic correctness is proved in a similar way as in the bivariate case (see Proposition 4.6), while also using Corollary 2.2. The sparse evaluation and interpolation at a geometric sequence require

operations. The univariate square-free factorizations can be done in time

, using [103] and [31, Theorem 14.23].

5.4.A pathological example

It is well known that a divisor of a sparse polynomial can have far more terms than the polynomial itself, because of the identity

and the possibility to replace

by any other sparse polynomial. For this reason, many problems on sparse polynomials turn out to be NP-hard [86, 87]. In a way that has been made precise in [25], this fundamental example is actually the main cause for such hardness results.

The example (5.1) is less dramatic if we consider sparse polynomials for which the total degree is much smaller than the total number of terms, which is often the case in practice. But even then, it still has some unpleasant consequences. Recall from [25] that

for coprime

with

. This gcd contains exactly

terms. Such gcds can also occur during content-free factorizations. For instance, the content of

Then

and

. The size of each individual factor in any irreducible factorization of

is bounded by

, which is sharp with respect to

. However, the content

has size

. This means that content-free factorization as a preparation step (before launching a “more expensive” factorization method) is a bad idea on this particular example.

6.Factoring using iterated Hensel lifting

Let

be a content-free polynomial and recall that any factor of

is again content-free. Let

be random non-zero elements of

. For any

and

, we define

In this section, we describe several algorithms for the factorization of

, following a similar recursive approach as for the computations of gcds in Section 3.1. This time, we start with the computation of a non-trivial factorization of the bivariate polynomial

. From this, we will recursively obtain non-trivial factorizations of

. This will in particular yield a non-trivial factorization of

. We will separately study the cases when

has a factorization into two or more coprime factors.

6.1.Faithful projections

In order to reconstruct factorizations of

from factorizations of

, it is important that the projection

be sufficiently generic. For

, let us denote by

the projection

. We say that

is faithful for

for

As usual, we assume that

are independently chosen at random from a subset of

of size at least

Proof. Given

, let

. Then

if and only if

, which happens with probability at most

Now consider a factorization

such that

are pairwise coprime. We say that

is faithful for this factorization if

is faithful for

and

are pairwise coprime.

Proof. This directly follows from Lemma 3.2 (while using that

) and the previous lemma (using induction on

While faithful mappings preserve coprime factorizations in a predictable way, it may still happen that an irreducible polynomial is projected to a reducible one. In fact this may even happen with probability

for a random choice of

Then

is irreducible, but

whenever

is a square in

. If

is a random element of

, then this happens with probability

. (Note that we may replace

by any other irreducible polynomial that involves both

and

This phenomenon is not so much of a problem if we want to recursively Hensel lift a coprime factorization

for which we know that there exist

with

and

(see Algorithms 6.1 and 6.2 below). However, any algorithm for Hensel lifting will fail if

is irreducible or, more generally, if no such

and

exist (Algorithm 6.3 below for the irreducible factorization of

may therefore fail on Example 6.3).

6.2.Lifting a coprime decomposition into two factors

Consider a non-trivial factorization

, where

are coprime. Assume that

is faithful for this factorization. Let us show how to recover

and

from

and

, for

With high probability, we have

, and the polynomials

and

are again coprime. It follows that each factorization

We now recover

and

from

and

, respectively, using sparse interpolation. In fact, if

, then we may interpolate

and recover

using one exact division.

Moreover, we may exploit the assumption that the supports of

and

with respect to

coincide with the supports of

and

. In the algorithm below, this explains why

evaluation points indeed suffice:

Input:	and coprime with
Output:	with , , and , whenever such a factorization exists (we raise an error if this is not the case) and is faithful for this factorization

Compute by recursively applying the algorithm to

Permute if necessary to ensure that

Compute , , and for , using sparse evaluation

Compute with using Hensel lifting (Section 4.3)

Raise an error if this fails

Recover from using sparse interpolation

Remark 6.4. The faithfulness assumption implies that

. If

is small, then this support bound is reasonably tight. Consequently, we may use sparse interpolation with known supports in order to recover

Proof. Assume that we correctly recovered

and

and let us investigate the probability that the results

and

are correct.

Let us first examine the probability that the Hensel lifting method from Section 4.3 works correctly. This is the case whenever the following three conditions are satisfied for

with

. We have

since

and

are coprime. The first condition fails for a given

if and only if

vanishes under the substitutions

. This happens with probability at most

for some

, by Corollary 2.2. Let us next consider

where

are formal indeterminates. Since

is content-free, we have

. Now

whenever

does not vanish under the substitutions

. The second condition therefore fails for some

with probability at most

, by Corollary 2.2 and using the fact that

. Finally, with probability at least

, we have

Let us now return to the correctness of the overall algorithm. Concerning the top-level call, we have the following:

Altogether, we correctly determine

and

with probability at least

. Since there are

recursive levels and since the sizes of the supports of

, and

are all smaller than those of

, and

, the probabilistic correctness claim follows.

As to the complexity bound, let us now study the cost when not counting the recursive calls:

We conclude that the cost is

, when not counting the recursive calls. Since there are

recursive levels and since the sizes of the supports of

, and

are all smaller than those of

, and

, the complexity bound follows.

Example 6.6. Let

, where

are coprime. Now consider the polynomial

. When factoring

using Algorithm 6.1, the last

lifting steps with respect to

all operate on a projection

of size

. This is an example when the overhead

in the complexity bound is sharp. However, the support of

is very special, since it is contained in an affine subspace of dimension

. It is easy to detect this situation and directly obtain the exponents in the last

variables as affine combinations of the known exponents in the first

variables.

More generally, it may happen that the fibers of the support under the projection with respect to the last

variables are all small. It would be interesting to develop specific optimizations for this situation, e.g. directly apply sparse interpolation on the fibers.

6.3.Simultaneous lifting of multiple factors

The lifting algorithm generalizes in a straightforward way to the case when

is a product of more than two factors. This time, we wish to recover a factorization

from

, assuming that

(hence

) are pairwise coprime.

Input:	, irreducible and pairwise coprime , , , and , such that
Output:	irreducible and pairwise coprime with and , whenever such a factorization exists (we raise an error if this is not the case) and is faithful for this factorization

Recursively apply the algorithm to to compute

Compute , for , using sparse evaluation

Compute with using Hensel lifting

Raise an error if this fails

Recover from the using sparse interpolation, for

Remark 6.7. More precisely, we compute

using binary powering and

from

using dense bivariate

-th root extraction, for

. As an optimization, we may first sort the factors of

such that

and then compute

using exact division.

Proof. The proof is similar to that of Algorithm 6.1. Assume now that we correctly recovered

and let us investigate the probability that the results

are correct. To apply the Hensel lifting method from Section 4.3, the following three conditions must be satisfied for

The analysis is similar as in the proof of Theorem 6.5, except that (6.1) now becomes

and the probability that

vanishes for one of the

substitutions is now bounded by

, since

. Using the fact there are

recursive levels, this completes the probabilistic correctness proof.

For the cost analysis, we again start by analyzing the cost of the algorithm without the recursive calls:

Summing these costs and multiplying with the recursive depths

yields the desired complexity bound.

Example 6.9. There are cases when the recursive application of Algorithm 6.1 may lead to intermediate expression swell. For instance, for

and

, consider the polynomials

Then a direct factorization of

into irreducibles can be done in time

, whereas recursively factoring out

and then

may lead us to consider intermediate expressions like

of size

Remark 6.10. The ideas behind Algorithm 6.2 can readily be adapted to

-adic lifting, in order to factor a polynomial

with rational coefficients. In that case, we start with a factorization

for some sufficiently large prime

. The goal is to lift this into a factorization of

for some sufficiently large

and then to recover a factorization of

using rational number reconstruction [31, Section 5.10]. The relevant analogues of

and

are

and

, where the prime number

plays a similar role as the formal indeterminate

. The bivariate Hensel lifting in

is replaced by traditional univariate Hensel lifting in

6.4.Irreducible factorization

In lucky cases, given an irreducible factorization

, random choices of

give rise with high probability to an irreducible factorization

. For such

, we may use the following algorithm to compute its irreducible factorization:

Input: a content-free polynomial

Output: an irreducible factorization or an exception

If , then return the result of a univariate irreducible factorization

Let be random elements of

Compute an irreducible factorization of

Apply Algorithm 6.2 to and this factorization

Remark 6.12. We have seen in Section 5.4 that the size of the content of a multivariate polynomial can be much larger than the size of the polynomial itself, in pathological cases. In such cases, we may wish to avoid content-free factorization as a first step of a general algorithm for irreducible factorization. This can be achieved by adapting Algorithm 6.2 so as to factor out the content

at the top-level and performing the recursive call and bivariate lifting to

instead of

. Incorporating the content-free factorization directly into the main algorithm in this way yields a similar complexity bound as in Theorem 6.11, but with

instead of

Remark 6.13. An alternative approach for reducing to the content-free case, which also works for the algorithms from Section 7 below, is to tag

according to a regularizing weight

. This has the advantage of “exposing” the entire factorization of

with respect to the tagging variable

. Of course, this requires to replace the degree

in the complexity bound by

Remark 6.14. The problem from Remark 6.3 can be partially remedied by amending the bivariate lifting step: instead of raising an error when

, we may continue the Hensel lifting to a higher precision (the

are now power series in

) and determine which factors need to be combined in order to obtain a factorization of

(using the techniques from [71]). Doing this for all

leads to a finest partition

such that

is a polynomial factor of

for

and

. We next apply the sparse interpolation to the polynomials

instead of the series

. This may require more than

interpolation points, so we also double

until the interpolation step succeeds. Combined with Remark 6.12, this yields to an algorithm of complexity

Unfortunately, we have not yet been able to perform a clean probabilistic analysis for this method. Indeed, we still need to prove that the algorithm terminates in reasonable expected time for irreducible polynomials. Now consider the example

Whenever

with

and

is a square in

, the projected polynomial

is not irreducible. In principle, this only happens with probability

, so this is unlikely to happen for all

values of

. However, we were unable so far to prove a general bound. We also note that this probabilistic argument is different from traditional Hilbert–Bertini style arguments [97, Section 6.1]; see also [71] and [72, Chapître 6].

Remark 6.15. Another way to turn Algorithm 6.3 into a method that always succeeds with high probability would be to apply random monomial transformation to

. For instance, after the change of variables

, we have

for the polynomial

from example (6.2). With high probability, Algorithm 6.3 returns the correct irreducible factorization of

after this rewriting.

7.Factoring using projective Hensel lifting

The iterative approach from Sections 3.1 and 6 has the disadvantage of introducing a dependence on the dimension

in the complexity bounds. In the case of gcd computations, the idea of regularizing weights allowed for the more direct approach in Section 3.2.

Unfortunately, direct adaptations of this approach to factorization only work in very special cases, because it is unclear how to recombine the factors found at different evaluation points. One notable exception is when

factors, say, into two irreducible factors of degrees one and two in one of the variables; in that case, we can recombine based on degree information.

In this section, we will develop another “direct” approach for the factorization of

that avoids iteration on the dimension

. Except for precomputations in the algorithm from Section 7.7, the algorithms in this section do not rely on regularizing weights, but exploits more subtle properties of the Newton polytope of

. We will present three versions of our approach in order to deal with Newton polytopes of increasing complexity. For simplicity, we only present them in the case when

is the product of two factors. As we shall see in Section 7.6, even the last and most elaborate version of our approach is not fully general. We will conclude with a theoretical variant that is less efficient in practice, but which has the advantage of providing a full algorithm for the computation of irreducible factorizations.

7.1.A favorable special case

where

are irreducible and distinct and both depend on the “main” variable, say,

. In particular, this implies that

is content-free and square-free. Assume also that the following conditions hold:

by requiring that

is monic. From now on we will always assume that we have done this.

Under these assumptions, we claim that

and

can be computed from

and

using Hensel lifting and sparse interpolation. In order to see this, let

be an admissible ratio or an FFT ratio. For each

, we define bivariate polynomials

Then

and

are coprime and we have

. This allows us to compute

and

from

, and

using Hensel lifting. For sufficiently large

, with

, we may then use sparse interpolation to recover

from

. This leads to the following algorithm:

Input:	such that there exists a factorization (7.1) that satisfies H1 and H2; we assume that and are given
Output:

using bivariate Hensel lifting for , ,

If yield through sparse interpolation, then return

Proof. The correctness is clear. Let us analyze the cost of the main steps of the algorithm:

Adding up these costs, the complexity bound follows. The algorithm is correct as long as the final division succeeds, which is the case with probability at least

7.2.Single slope bivariate Hensel lifting

The assumptions H1 and H2 are fairly strong: they are never satisfied both if

is a homogeneous polynomial. Now the assumptions H1 and H2 essentially allow us to perform the bivariate Hensel lifting using the method from Section 4.3. The problem that we face in order to extend the approach from Algorithm 7.1 is that random shifts

are not allowed for the bivariate Hensel lifting from Section 4.3: we only have an efficient algorithm to compute

and

are given, not

and

Instead, we need a way to generalize the algorithm from Section 4.3 to the case when the Newton polygon of

is non-trivial. In this subsection, we start with the simplest “single slope” case. Since this material is a supplement to Section 4.3, we adopt the notation from there.

Assume that

has a non trivial factorization

with

. We define the Newton polygon of

to be the “lower border” of its Newton polytope:

It is the union of a finite number of edges

between

with

, and

. For each edge

, there is a corresponding weight

such that

Assume that

has a single edge, let

, and let

and

be coprime such that

. Now consider

Then

and

. Moreover,

, and we have similar expressions for

and

. If

and

are known and if their transforms

and

are coprime, then this enables us to compute

and

by applying the Hensel lifting method from Section 4.3 to

, and

7.3.Single slope factorization

Let us now return to the factorization problem from Section 7.1 and let us show how to generalize Algorithm 7.1 using the material from the previous subsection.

Let

and let

be the lower boundary of its convex hull. With high probability, we have

and

, for any given

. The last edge of

joins the points

and

for

with

. Let

with

and

, where

, and

. Now consider the assumptions

The hypotheses H1 and H2 correspond to the special case when

. If

and

are known (e.g. through the recursive factorization of

is not

-homogeneous), then the algorithm from Section 7.2 allows us to Hensel lift the factorization

into a factorization

. This leads to the following generalization of Algorithm 7.1:

Input:	such that there exists a factorization (7.1) that satisfies H1' and H2', together with and such that
Output:

using bivariate Hensel lifting from Section 7.2 for , ,

If yield through sparse interpolation, then return

Proof. The proof is similar to the proof of Theorem 7.1. The main change is that the generalized algorithm also requires the evaluations of

and

on our sequence. This takes

additional operations, which is absorbed by the complexity bound.

Note that the assumptions H1' and H2' are actually slightly more liberal than the statement that the Newton polygon should have a single edge. For this reason, Algorithm 7.2 is very likely to apply as soon as there exists a variable

for which

is small (of course, we may then assume

modulo a permutation of variables). This in particular the case when

, unless

is not square-free (one may need to replace

and multiply

with

). Note that a single “good” variable

suffices.

Remark 7.3. Both Algorithms 6.2 and 7.2 involve recursive factorizations of polynomials in less variables. However, the supports of the recursive calls are obtained through projection for Algorithm 6.2 and through restriction for Algorithm 7.2. The second case is often more favorable in the sense that the recursive supports decrease faster in size.

Remark 7.4. Consider the case when

and

. Then

coincides with

, so the requested factorization of

is not simpler than the desired factorization of

. This is due to the fact that the exponents of

all lie in a linear subspace of

. In particular, the degree in

of any term of

can be determined as a function of the degrees of the other variables. In particular, we may recover a factorization of

from a factorization of

. Using a straightforward induction on

, this allows us to ensure that

has strictly less terms than

7.4.Multiple slope bivariate Hensel lifting

As soon as

, it can happen that every factor of

has a Newton polygon with at least two edges. In order to deal with this situation, the algorithm from Section 7.2 may be further generalized to accommodate Newton polygons with multiple slopes. In order to explain how, we again adopt the notation from Sections 4.3 and 7.2.

So assume that

and

, that

has an arbitrary number of edges, and that

and

are known for each edge

. Assume also that

and

are coprime for

, where

is the normalization mapping with

. We will say that

is a coprime edge factorization of

. The aim of this subsection is to lift these factorizations efficiently into the global factorization

It is well known that a variant of the Newton polygon method can be used in order to factor

over the ring

of Laurent series. An efficient algorithm for this task was described in [52]. One important subtask is distinct slope factorization: each edge

gives rise to a unique monic factor

such that

for

, and the natural generalization of the notation

. The methods from [52] allow for the efficient computation of the

: it takes time

to compute

at order

for

Now for each

, the known factor

induces a factor of

that can be Hensel lifted into the factor

, by adapting the techniques from Section 4.3 to Laurent series coefficients. For some

, we then have

. We may determine

in a similar way as in Section 4.3. The other factor

can be obtained using one bivariate division. The total cost of this method to compute

and

is bounded by

7.5.Multiple slope factorization

We are now in a position to generalize Algorithm 7.2 to the case when the bivariate Hensel lifting is done using the multiple slope method from the previous subsection. The general algorithm is fairly technical, so we found it more instructive to illustrate the main ideas on an example. Let

whence

is the closed triangle with vertices

, and

. Its lower boundary consists of the edge from

and the edge from

, which correspond to the weights

and

respectively. The generalized algorithm starts with the recursive factorizations of

and

, which yields

However, these factorizations do not tell us from which factors of

they originate: the factorization of

could have been of the form

, with

and

. In order to obtain the correct matches, the next idea is to consider the bivariate factorization of

for some sufficiently generic values of

and

. For instance, for

and

, we obtain

From these computations, we deduce that the factor

comes from the same factor of

as the factor

. We say that we have managed to match the factors

and

of the different slopes. In other words, if

is a non-trivial factor of

, then, up to constant scaling, we now know that either

For any

and

, let

. Having completed the above matching, we may now use the lifting algorithm from Section 7.4 to deduce

and

from

, and

. Doing this for sufficiently many

, we may finally recover

and

using sparse interpolation.

7.6.A torture example

The approach from Section 7.5 is fairly general. However, as we will show now, it is possible to construct pathological examples that can still not be treated with this method.

Example 7.5. Consider a polynomial

and a monomial

such that

lies in the interior of the Newton polytope

. For instance, we may take

for pairwise distinct

. This polynomial cannot be factored using the techniques from this section, so far.

7.7.Factor matching through homotopies

In fact, it is not so much the factorization of the bivariate

polynomials that is a problem: instead of Hensel lifting, we may very well rely on Theorem 4.7 (this only increases the exponent in

in Theorems 7.1 and 7.2, which is subdominant anyway if

). The true problem is matching corresponding factors of

for different

, especially in the most difficult cases like Example 7.5. We conclude this section with an approach that completely solves this problem, modulo a polynomial overhead in

is irreducible, then

is irreducible with high probability, by Theorem 1.1, and so are

and

. For general

, this also means that the factors of

, and

are in effective one to one correspondence, with high probability. Moreover, by construction,

Using this relation, we may therefore match corresponding factors of

and

with high probability.

Having solved the matching problem, we still need to ensure that the factorizations of the

are normalized in a suitable way such that we can recover the factors of

using sparse interpolation. For this, let

be a regularizing weight for

with

. Let

be yet another indeterminate and consider

By construction, any factor

has leading coefficient

with respect to

, for some

, and some monomial

. Thus, any such factor

can be normalized by dividing out the constant

. Now consider the unique factorization of

into a product of normalized factors times a non-zero scalar in

. This yields a factorization of

by setting

The factorization of

can be lifted to the factorization of

, which is, in turn, lifted to the factorization of

, and so on. An irreducible factorization of

is recovered via interpolation from the factorizations of

for sufficiently large

. Finally, to obtain a factorization of

from a factorization of

, we can set

, and apply the map

Output: an irreducible factorization of

Let be a regularizing weight for

Compute an irreducible factorization , normalized as above

into an irreducible factorization

Deduce from via root extraction

Try to determine from () using sparse interpolation

If successful, then set , and return

Proof. We factor

as a dense polynomial in four variables of degree

(after division by a suitable power of

) using the algorithm from [69, Proposition 5]. This requires

operations in

and the probability of success is at least

Let

be an irreducible factorization of

. By our assumption on the characteristic, the factors

are all separable. We will apply [72, Théorème 6.9] for each of the points in our geometric progression. By using Corollary 2.2 instead of the Schwartz-Zippel lemma in the proof of [72, Théorème 6.9], we deduce that the specializations

, and thus

are irreducible for

and

, with probability at least

. Under this assumption, and modulo reordering the factors, the computed

are of the form

for suitable scaling factors

that do not depend on

. The check whether the computed factorization of

is correct reports a false positive with probability at most

, by Remark 2.4 or the Schwarz-Zippel lemma.

Let us now analyze the complexity. We first observe that

, since

when

consists of a single monomial of total degree

. Using [2, Theorem 5] in a recursive manner, we may compute

from

in time

. We saw above that the factorization of

requires at most

operations in

. The computation of the specializations

for

requires

further operations. The Hensel lifting of the

can be done in time

using Proposition 4.4 and evaluation-interpolation in the remaining variable. The

Hensel lifts succeed with probability at least

. Adding up, we obtain the desired complexity bound, as well as the claimed bound for the probability of success.

Remark 7.7. Due to the

overhead, Algorithm 7.3 is mainly of theoretical interest. It is plausible that (7.2) can be replaced by a less, but still sufficiently, generic formula. This would reduce the overhead in

. In practice, one may use Algorithm 7.3 as a last resort, in the unlikely case that all other strategies from Sections 6 and 7 fail.

Bibliography

[1]: F. Abu Salem, S. Gao, and A. G. B. Lauder. Factoring polynomials via polytopes. In Proc. ISSAC '04, pages 4–11. New York, NY, USA, 2004. ACM.
[2]: A. V. Aho, K. Steiglitz, and J. D. Ullman. Evaluating polynomials on a fixed set of points. SIAM Journ. of Comp., 4:533–539, 1975.
[3]: M. Alberich-Carramiñana, J. Guàrdia, E. Nart, A. Poteaux, J. Roé, and M. Weimann. Polynomial factorization over henselian fields. Technical Report, arXiv, 2022. https://doi.org/10.48550/arXiv.2207.02139.
[4]: K. Belabas, M. van Hoeij, J. Klüners, and A. Steel. Factoring polynomials over global fields. J. Théor. Nombres Bordeaux, 21(1):15–39, 2009.
[5]: M. Ben-Or and P. Tiwari. A deterministic algorithm for sparse multivariate polynomial interpolation. In Proc. ACM STOC '88, pages 301–309. New York, NY, USA, 1988.
[6]: E. R. Berlekamp. Factoring polynomials over finite fields. Bell System Tech. J., 46:1853–1859, 1967.
[7]: E. R. Berlekamp. Factoring polynomials over large finite fields. Math. Comp., 24:713–735, 1970.
[8]: L. Bernardin. Factorization of multivariate polynomials over finite fields. PhD thesis, ETH, Zürich, 1999.
[9]: D. J. Bernstein. Scaled remainder trees. Available from https://cr.yp.to/arith/scaledmod-20040820.pdf, 2004.
[10]: E. Bertini. Introduzione alla geometria proiettiva degli iperspazi con appendice sulle curve algebriche e loro singolarità. Giuseppe Principato, 1923.
[11]: A. Borodin and R. T. Moenck. Fast modular transforms. Journal of Computer and System Sciences, 8:366–386, 1974.
[12]: A. Bostan, G. Lecerf, and É. Schost. Tellegen's principle into practice. In Proc. ISSAC '03, pages 37–44. Philadelphia, USA, August 2003.
[13]: A. Bostan and É. Schost. Polynomial evaluation and interpolation on special sets of points. J. of Complexity, 21(4):420–446, August 2005. Festschrift for the 70th Birthday of Arnold Schönhage.
[14]: P. Bürgisser, M. Clausen, and M. A. Shokrollahi. Algebraic complexity theory. Springer-Verlag, Berlin, 1997.
[15]: D. G. Cantor and E. Kaltofen. On fast multiplication of polynomials over arbitrary algebras. Acta Informatica, 28:693–701, 1991.
[16]: D. G. Cantor and H. Zassenhaus. A new algorithm for factoring polynomials over finite fields. Math. Comp., 36(154):587–592, 1981.
[17]: T. Chen and M. Monagan. The complexity and parallel implementation of two sparse multivariate Hensel lifting algorithms for polynomial factorization. In Proc. CASC 2020, pages 150–169. Springer-Verlag, 2020.
[18]: T. Chen and M. Monagan. Factoring multivariate polynomials represented by black boxes: a maple + C implementation. JSC, 16, 2022.
[19]: G. Chèze and G. Lecerf. Lifting and recombination techniques for absolute factorization. J. of Complexity, 23(3):380–420, 2007.
[20]: A. L. Chistov. Efficient factoring polynomials over local fields and its applications. In Proceedings of the International Congress of Mathematicians, Kyoto, Japan, 1990, volume 1, pages 1509–1519. Springer-Verlag, 1991.
[21]: G. E. Collins. Subresultants and reduced polynomial remainder sequences. JACM, 14(1):128–142, 1967.
[22]: A. Cuyt and W.-S. Lee. Sparse interpolation of multivariate rational functions. TCS, 412:1445–1456, 2011.
[23]: W. Decker, G. Lecerf, and G. Pfister. Absolute factorization for characteristic . Singular library absfact.lib, 2005.
[24]: J. Della Dora, C. Dicrescenzo, and D. Duval. A new method for computing in algebraic number fields. In Eurocal'85 (2), volume 174 of Lect. Notes in Comp. Science, pages 321–326. Springer, 1985.
[25]: M. Filaseta, A. Granville, and A. Schinzel. Irreducibility and greatest common divisor algorithms for sparse polynomials. London Math. Soc. Lect. Note Series, 352:155–176, 2008.
[26]: D. Ford. On the computation of the maximal order in a Dedekind domain. PhD thesis, Ohio State University, 1978.
[27]: A. Fröhlich and J. C. Shepherdson. On the factorisation of polynomials in a finite number of steps. Math. Z., 62:331–334, 1955.
[28]: A. Fröhlich and J. C. Shepherdson. Effective procedures in field theory. Philos. Trans. Soc. London. Ser. A., 248:407–432, 1956.
[29]: S. Gao. Factoring multivariate polynomials via partial differential equations. Math. Comp., 72(242):801–822, 2003.
[30]: J. von zur Gathen. Irreducibility of multivariate polynomials. J. Comp. System Sci., 31:225–264, 1985.
[31]: J. von zur Gathen and J. Gerhard. Modern Computer Algebra. Cambridge University Press, New York, NY, USA, 3rd edition, 2013.
[32]: J. von zur Gathen and E. Kaltofen. Factoring sparse multivariate polynomials. J. Comput. System Sci., 31(2):265–287, 1983. Proc. FOCS.
[33]: J. von zur Gathen and E. Kaltofen. Factorization of multivariate polynomials over finite fields. Math. Comp., 45(171):251–261, 1985.
[34]: P. Gianni and B. Trager. Square-free algorithms in positive characteristic. AAECC, 7(1):1–14, 1996.
[35]: M. Giesbrecht and D. S. Roche. Diversification improves interpolation. In Proc. ISSAC '11, pages 123–130. New York, NY, USA, 2011. ACM.
[36]: B. Grenet. Bounded-degree factors of lacunary multivariate polynomials. JCS, 75:171–192, 2016.
[37]: B. Grenet, J. van der Hoeven, and G. Lecerf. Randomized root finding over finite fields using tangent Graeffe transforms. In Proc. ISSAC '15, pages 197–204. New York, NY, USA, 2015. ACM.
[38]: J. Guàrdia, J. Montes, and E. Nart. Newton polygons of higher order in algebraic number theory. Trans. Amer. Math. Soc., 364(1):361–416, 2012.
[39]: Jr. H. W. Lenstra. Finding small degree factors of lacunary polynomials. In In Number theory in progress, volume 1, pages 267–276. De Gruyter, Zakopane-Kościelisko, 1999.
[40]: D. Harvey and J. van der Hoeven. Polynomial multiplication over finite fields in time . Technical Report, HAL, 2019. http://hal.archives-ouvertes.fr/hal-02070816.
[41]: K. Hensel. Neue Grundlagen der Arithmetik. J. Reine Angew. Math., 127:51–84, 1904.
[42]: D. Hilbert. Über die Irreducibilität ganzer rationaler Functionen mit ganzzahligen Coefficienten. J. reine angew. Math., 110:104–129, 1892.
[43]: M. van Hoeij. Factoring polynomials and the knapsack problem. Journal of Number theory, 95(2):167–189, 2002.
[44]: J. van der Hoeven. Faster Chinese remaindering. Technical Report, HAL, 2016. https://hal.archives-ouvertes.fr/hal-01403810.
[45]: J. van der Hoeven. Probably faster multiplication of sparse polynomials. Technical Report, HAL, 2020. https://hal.archives-ouvertes.fr/hal-02473830.
[46]: J. van der Hoeven. The Jolly Writer. Your Guide to GNU TeXmacs. Scypress, 2020.
[47]: J. van der Hoeven and G. Lecerf. On the bit-complexity of sparse polynomial multiplication. JSC, 50:227–254, 2013.
[48]: J. van der Hoeven and G. Lecerf. Sparse polynomial interpolation. Exploring fast heuristic algorithms over finite fields. Technical Report, HAL, 2019. https://hal.science/hal-02382117v1.
[49]: J. van der Hoeven and G. Lecerf. Directed evaluation. J. of Complexity, 60, 2020. Article ID 101498.
[50]: J. van der Hoeven and G. Lecerf. Fast multivariate multi-point evaluation revisited. J. of Complexity, 56, 2020. Article ID 101405, 38 pages.
[51]: J. van der Hoeven and G. Lecerf. On sparse interpolation of rational functions and gcds. ACM SIGSAM Commun. Comput. Algebra, 55(1):1–12, 2021.
[52]: J. van der Hoeven and G. Lecerf. Approximate contact factorization of germs of plane curves. Technical Report, HAL, 2022. https://hal.archives-ouvertes.fr/hal-03745581.
[53]: J. van der Hoeven and G. Lecerf. Univariate polynomial factorization over large finite fields. AAECC, 2022. https://doi.org/10.1007/s00200-021-00536-1.
[54]: J. Hu and M. B. Monagan. A fast parallel sparse polynomial GCD algorithm. In S. A. Abramov, E. V. Zima, and X.-S. Gao, editors, Proc. ISSAC '16, pages 271–278. ACM, 2016.
[55]: Q.-L. Huang and X.-S. Gao. Bit complexity of polynomial gcd on sparse representation. Preprint available from https://arxiv.org/abs/2207.13874, 2022.
[56]: M. Javadi and M. Monagan. A sparse modular gcd algorithm for polynomial gcd computation over algebraic function fields. In Proc. ISSAC '07, pages 187–194. 2007.
[57]: E. Kaltofen. A polynomial reduction from multivariate to bivariate integral polynomial factorization. In Proc. STOC '82, pages 261–266. ACM, 1982.
[58]: E. Kaltofen. Polynomial-time reductions from multivariate to bi- and univariate integral polynomial factorization. SIAM J. Comput., 14(2):469–489, 1985.
[59]: E. Kaltofen. Sparse Hensel lifting. In Proc. EUROCAL '85, pages 4–17. Berlin, Heidelberg, 1985. Springer.
[60]: E. Kaltofen. Factorization of polynomials given by straight-line programs. In S. Micali, editor, Randomness and Computation, volume 5 of Advances in Computing Research, pages 375–412. JAI Press Inc., Greenwhich, Connecticut, 1989.
[61]: E. Kaltofen. Effective noether irreducibility forms and applications. J. Computer and System Sciences, 50(2):274–295, 1995.
[62]: E. Kaltofen and G. Lecerf. Handbook of finite fields, chapter Factorization of multivariate polynomials, pages 383–392. Discrete Mathematics and its Applications. Chapman and Hall/CRC, 2013.
[63]: E. Kaltofen and M. B. Monagan. On the genericity of the modular polynomial gcd algorithm. In Proc. ISSAC '99, pages 59–66. 1999.
[64]: E. Kaltofen and V. Shoup. Subquadratic-time factoring of polynomials over finite fields. Math. Comp., 67(223):1179–1197, 1998.
[65]: E. Kaltofen and B. M. Trager. Computing with polynomials given by black boxes for their evaluations: greatest common divisors, factorization, separation of numerators and denominators. JSC, 9(3):301–320, 1990.
[66]: E. Kaltofen and L. Yagati. Improved sparse multivariate polynomial interpolation algorithms. In ISSAC '88: Proceedings of the International Symposium on Symbolic and Algebraic Computation, pages 467–474. Springer Verlag, 1988.
[67]: K. S. Kedlaya and C. Umans. Fast polynomial factorization and modular composition. SIAM J. Comput., 40(6):1767–1802, 2011.
[68]: G. Lecerf. Sharp precision in Hensel lifting for bivariate polynomial factorization. Mathematics of Computation, 75:921–933, 2006.
[69]: G. Lecerf. Improved dense multivariate polynomial factorization algorithms. JSC, 42(4):477–494, 2007.
[70]: G. Lecerf. Fast separable factorization and applications. Applicable Algebra in Engineering, Communication and Computing, 19(2), 2008.
[71]: G. Lecerf. New recombination algorithms for bivariate polynomial factorization based on Hensel lifting. AAECC, 21(2):151–176, 2010.
[72]: G. Lecerf. Journées Nationales de Calcul Formel (2013), volume 3 of Les cours du CIRM, chapter Factorisation des polynômes à plusieurs variables. CEDRAM, 2013. Exp. No. 2, 85 pages, https://ccirm.centre-mersenne.org/articles/10.5802/ccirm.18.
[73]: G. Lecerf. On the complexity of the Lickteig–Roy subresultant algorithm. JSC, 92:243–268, 2019.
[74]: A. K. Lenstra, H. W. Lenstra, and L. Lovász. Factoring polynomials with rational coefficients. Math. Ann., 261:515–534, 1982.
[75]: M. Monagan and B. Tuncer. Using sparse interpolation in Hensel lifting. In Proceedings of CASC 2016, volume 9890 of LNCS, pages 381–400. Springer, 2016.
[76]: M. Monagan and B. Tuncer. Factoring multivariate polynomials with many factors and huge coefficients. In Computer Algebra in Scientific Computing, pages 319–334. Cham, 2018. Springer International Publishing.
[77]: M. Monagan and B. Tuncer. Sparse multivariate Hensel lifting: a high-performance design and implementation. In Proc. ICMS 2018, volume 10931 of LNCS, pages 359–368. Springer, 2018.
[78]: M. Monagan and B. Tuncer. The complexity of sparse Hensel lifting and sparse polynomial factorization. Journal of Symbolic Computation, 99:189–230, 2020.
[79]: J. Montes. Polígonos de Newton de orden superior y aplicaciones aritméticas. PhD thesis, Universitat de Barcelona, Spain, 1999.
[80]: G. Moroz. New data structure for univariate polynomial approximation and applications to root isolation, numerical multipoint evaluation, and other problems. In Proc. IEEE FOCS '21. Denver, United States, 2022.
[81]: D. R. Musser. Multivariate polynomial factorization. JACM, 22(2):291–308, 1975.
[82]: A. Novocin. Factoring Univariate Polynomials over the Rationals. PhD thesis, Florida State University at Tallahassee, 2008.
[83]: A. Novocin, D. Stehlé, and G. Villard. An LLL-Reduction algorithm with quasi-linear time complexity: extended abstract. In Proc. ACM STOC '11, pages 403–412. New York, NY, USA, 2011.
[84]: V. Y. Pan. Univariate polynomials: nearly optimal algorithms for numerical factorization and root-finding. JSC, 33(5):701–733, 2002.
[85]: C. H. Papadimitriou. Computational Complexity. Addison-Wesley, 1994.
[86]: D. A. Plaisted. Sparse complex polynomials and polynomial reducibility. J. of Computer and System Sciences, 14(2):210–221, 1977.
[87]: D. A. Plaisted. New NP-hard and NP-complete polynomial and integer divisibility problems. TCS, 31(1-2):125–138, 1984.
[88]: R. Prony. Essai expérimental et analytique sur les lois de la dilatabilité des fluides élastiques et sur celles de la force expansive de la vapeur de l'eau et de la vapeur de l'alkool, à différentes températures. J. de l'École Polytechnique Floréal et Plairial, an III, 1:24–76, 1795. Cahier 22.
[89]: D. S. Roche. What can (and can't) we do with sparse polynomials? In Proc. ISSAC '18, pages 25–30. New York, NY, USA, 2018. ACM.
[90]: W. M. Ruppert. Reduzibilität ebener kurven. J. Reine Angew. Math., 396:167–191, 1986.
[91]: K. Ryan and N. Heninger. Fast practical lattice reduction through iterated compression. In Advances in Cryptology – CRYPTO 2023, pages 3–36. 2023.
[92]: T. Sasaki and M. Sasaki. A unified method for multivariate polynomial factorizations. Japan J. Indust. Appl. Math., 10(1):21–39, 1993.
[93]: N. Saxena. Closure of algebraic classes under factoring. https://www.cse.iitk.ac.in/users/nitin/talks/Sep2023-paris.pdf, 2023. Talk at Recent Trends in Computer Algebra, Institut Henri Poincar´e, Paris.
[94]: Th. Schönemann. Von denjenigen Moduln, welche Potenzen von Primzahlen sind. J. Reine Angew. Math., 32:93–105, 1846. See §59.
[95]: A. Schönhage. Asymptotically fast algorithms for the numerical multiplication and division of polynomials with complex coefficients. In J. Calmet, editor, EUROCAM '82: European Computer Algebra Conference, volume 144 of Lect. Notes Comp. Sci., pages 3–15. Marseille, France, April 1982. Springer.
[96]: J. T. Schwartz. Fast probabilistic algorithms for verification of polynomial identities. JACM, 27(4):701–717, 1980.
[97]: I. R. Shafarevich. Basic Algebraic Geometry 1. Springer, 3rd edition, 2013.
[98]: V. Shoup. On the deterministic complexity of factoring polynomials over finite fields. Information Processing Letters, pages 261–267, 1990.
[99]: A. Steel. Conquering inseparability: primary decomposition and multivariate factorization over algebraic function fields of positive characteristic. JSC, 40(3):1053–1075, 2005.
[100]: P. S. Wang. An improved multivariate polynomial factoring algorithm. Math. Comp., 32(144):1215–1231, 1978.
[101]: P. S. Wang and L. P. Rothschild. Factoring multivariate polynomials over the integers. Math. Comp., 29:935–950, 1975.
[102]: W. Wu, J. Chen, and Y. Feng. Sparse bivariate polynomial factorization. Science China Mathematics, 57:2123–2142, 2014.
[103]: D. Y. Y. Yun. On square-free decomposition algorithms. In Proc. SYMSAC '76, pages 26–35. New York, NY, USA, 1976. ACM.
[104]: H. Zassenhaus. On Hensel factorization, I. Journal of Number Theory, 1(3):291–311, 1969.
[105]: R. Zippel. Probabilistic algorithms for sparse polynomials. In Proc. EUROSAM '79, pages 216–226. Springer, 1979.
[106]: R. Zippel. Newton's iteration and the sparse Hensel algorithm (extended abstract). In Proc. SYMSAC '81, pages 68–72. New York, NY, USA, 1981. ACM.

1.Introduction

1.1.Motivation and main goals

1.2.Overview of univariate and bivariate factorization methods

1.3.Overview of multivariate factorization methods

1.4.Outline of our contributions

1.5.Notation

2.Preliminaries

2.1.Basic complexities

2.2.The Schwarz–Zippel lemma

2.3.Sparse polynomial interpolation

2.4.Sparse evaluation at geometric progressions

2.5.Newton polytopes

2.6.Laurent polynomials

2.7.Tagging

3.Multivariate gcd computations

3.1.Iterative computation of gcds

3.2.Gcd computations through regularizing weights

4.Bivariate factorization

4.1.Content-free factorization

4.2.Root extraction

4.3.Hensel lifting

4.4.Square-free factorization

4.5.General bivariate factorization

5.Efficient reductions

5.1.Content-free factorization

5.2.Root extraction

5.3.Square-free factorization

5.4.A pathological example

6.Factoring using iterated Hensel lifting

6.1.Faithful projections

6.2.Lifting a coprime decomposition into two factors

6.3.Simultaneous lifting of multiple factors

6.4.Irreducible factorization

7.Factoring using projective Hensel lifting

7.1.A favorable special case

7.2.Single slope bivariate Hensel lifting

7.3.Single slope factorization

7.4.Multiple slope bivariate Hensel lifting

7.5.Multiple slope factorization

7.6.A torture example

7.7.Factor matching through homotopies

Bibliography