We consider the problem of functional decomposition of univariate sparse polynomials. Suppose \(f\in {\mathsf {F}}[x]\) is a polynomial over a field \({\mathsf {F}}\) which has a non-trivial decomposition into \(g,h\in {\mathsf {F}}[x]\) such that f(x) = g(h(x)). Given f and \(\deg g\), our algorithm produces both decomposition factors \(g,h\in {\mathsf {F}}[x]\) in polynomial time in the number of nonzero terms in f, h, the degree of g, and \(\log \deg f\). Work by Zannier implies that \(\deg g\) is typically bounded quadratically by the number of non-zero terms of f, so that we can say the algorithm runs in polynomial time in the sparse sizes of the inputs and outputs. Our algorithm works when the characteristic of \({\mathsf {F}}\) is 0 or does not divide \(\deg g\). We also develop a heuristic and probabilistic classifier, based on random sampling of evaluations in a finite field, which can (on average) distinguish whether a given polynomial is decomposable or not, in sub-linear time in \(\deg f\). Finally, we give a generalization to sparse systems of linear differential equations.

We propose a randomized algorithm for multiplication in the ring of non-commutative polynomials Κ [x1,…,xn]{#948;1,…,δn}, where δ i=xi∂over∂ xi, dedicated to sparse inputs. The complexity of our algorithm is polynomial in the input size and on an a priori sparsity bound for the output.

We present a full analysis of the bit complexity of an efficient algorithm for the computation of at least one point in each connected component of a smooth real hypersurface. This is a basic and important operation in semi-algebraic geometry: it gives an upper bound on the number of connected components of a real hypersurface, and is also used in many higher level algorithms.

Consider the skew polynomial ring L[x; σ], where L is a field and σ is an automorphism of L of order r. We present two randomized algorithms for the multiplication of sparse skew polynomials in L[x;σ].

We present a new algorithm for the computation of the integer-linear decomposition of a multivariate polynomial. Such a decomposition is used in Ore-Sato theory and discrete creative telescoping, for example to detect applicability of Zeilberger's algorithm to a hypergeometric term. Our algorithm is quite straightforward, requiring only basic polynomial arithmetic along with efficient rational root finding. We present complete complexity analyses for both our and previous algorithms in the case of bivariate integer polynomials, and show that our method has a better theoretical performance. We also provide a Maple implementation which shows that our method is faster in practice than previous algorithms.

For any finite Galois field extension K/F, with Galois group G = Gal(K/F), there exists an element K whose orbit G · forms an F-basis of K. Such an is called a normal element and G · is a normal basis. We introduce a probabilistic algorithm for finding a normal element when G is either a finite abelian or a metacyclic group. The algorithm is based on the fact that deciding whether a random element K is normal can be reduced to deciding whether () K[G] is invertible (where ranges over all of G). Our algorithm requires a quadratic number of operations in the size of G for metacyclic G, and a slightly subquadratic number of operations for abelian G.

We consider the problem of computing the nearest matrix polynomial with a non-trivial Smith Normal Form. We show that computing the Smith form of a matrix polynomial is amenable to numeric computation as an optimization problem. Furthermore, we describe an effective optimization technique to find a nearby matrix polynomial with a non-trivial Smith form. The results are later generalized to include the computation of a matrix polynomial having a maximum specified number of ones in the Smith Form (i.e., with a maximum specified McCoy rank). We discuss the geometry and existence of solutions and how our results can used for a backwards error analysis. We develop an optimization-based approach and demonstrate an iterative numerical method for computing a nearby matrix polynomial with the desired spectral properties. We also describe the implementation of our algorithms and demonstrate the robustness with examples in Maple.

Matrix polynomials appear in many areas of computational algebra, control systems theory, differential equations, and mechanics, typically with real or complex coefficients. Because of numerical error and instability, a matrix polynomial may appear of considerably higher rank (generically full rank), while being very close to a rank-deficient matrix. "Close" is defined naturally under the Frobenius norm on the underlying coefficient matrices of the matrix polynomial. In this paper we consider the problem of finding the nearest rank-deficient matrix polynomial to an input matrix polynomial, that is, the nearest square matrix polynomial which is algebraically singular. We prove that such singular matrices at minimal distance always exist (and we are never in the awkward situation having an infimum but no actual matrix polynomial at minimal distance). We also show that singular matrices at minimal distance are all isolated, and are surrounded by a basin of attraction of non-minimal solutions. Finally, we present an iterative algorithm which, on given input sufficiently close to a rank-deficient matrix, produces that matrix. The algorithm is efficient and is proven to converge quadratically given a sufficiently good starting point. An implementation demonstrates the effectiveness and numerical robustness in practice.

We present a new Monte Carlo algorithm for the interpolation of a straight-line program as a sparse polynomial f over an arbitrary finite field of size q. We assume a priori bounds D and T are given on the degree and number of terms of f. The approach presented in this paper is a hybrid of the diversified and recursive interpolation algorithms, the two previous fastest known probabilistic methods for this problem. By making effective use of the information contained in the coefficients themselves, this new algorithm improves on the bit complexity of previous methods by a "soft-Oh" factor of T, log D, or log q.

In this paper, we present a new algorithm and an experimental implementation for factoring elements in the polynomial nth Weyl algebra, the polynomial nth shift algebra, and Zn-graded polynomials in the nth q-Weyl algebra.

We present algorithms to compute the Smith Normal Form of matrices over two families of local rings. The algorithms use the black-box model which is suitable for sparse and structured matrices. The algorithms depend on a number of tools, such as matrix rank computation over finite fields, for which the best-known time- and memory-efficient algorithms are probabilistic.

We consider the problem of interpolating an unknown multivariate polynomial with coefficients taken from a finite field or as numerical approximations of complex numbers. Building on the recent work of Garg and Schost, we improve on the best-known algorithm for interpolation over large finite fields by presenting a Las Vegas randomized algorithm that uses fewer black box evaluations. Using related techniques, we also address numerical interpolation of sparse polynomials with complex coefficients, and provide the first provably stable algorithm (in the sense of relative error) for this problem, at the cost of modestly more evaluations. A key new technique is a randomization which makes all coefficients of the unknown polynomial distinguishable, producing what we call a diverse polynomial. Another departure from most previous approaches is that our algorithms do not rely on root finding as a subroutine. We show how these improvements affect the practical performance with trial implementations.

The functional decomposition of polynomials has been a topic of great interest and importance in pure and computer algebra and their applications. The structure of compositions of (suitably normalized) polynomials f = g o h in Fq[x] is well understood in many cases, but quite poorly when the degrees of both components are divisible by the characteristic p. This work investigates the decomposition of polynomials whose degree is a power of p. An (equal-degree) i-collision is a set of i distinct pairs (g, h) of polynomials, all with the same composition and deg g the same for all (g, h). Abhyankar (1997) introduced the projective polynomials xn+ ax + b, where n is of the form (rm -- 1)/(r -- 1) and r is a power of p. Our first tool is a bijective correspondence between i-collisions of certain additive trinomials, projective polynomials with i roots, and linear spaces with i Frobenius-invariant lines.

We consider the problem of determining whether a t-sparse or lacunary polynomial f is a perfect power, that is, f=hr for some other polynomial h and positive integer r, and of finding h and r should they exist. We show how to determine if f is a perfect power in time polynomial in the size of the lacunary representation. The algorithm works over GF(q)[x] (at least for large characteristic) and over Z[x], where the cost is also polynomial in the log of the infinity norm of f. Subject to a conjecture, we show how to find h if it exists via a kind of sparse Newton iteration, again in time polynomial in the size of the sparse representation. Finally, we demonstrate an implementation using the C++ library NTL.

Efficient block projections of non-singular matrices have recently been used by the authors in [10] to obtain an efficient algorithm to find rational solutions for sparse systems of linear equations. In particular a bound ofO~(n2.5) machine operations is presented for this computation assuming that the input matrix can be multiplied by a vector with constant-sized entries using O~(n) machine operations. Somewhat more general bounds for black-box matrix computations are also derived. Unfortunately, the correctness of this algorithm depends on the existence of efficient block projections of non-singular matrices, and this was only conjectured.

We propose a new algorithm to find a rational solution to a sparse system of linear equations over the integers. This algorithm is based on a p-adic lifting technique combined with the use of block matrices with structured blocks. It achieves a sub-cubic complexity in terms of machine operations subject to a conjecture on the effectiveness of certain sparse projections. A LinBox-based implementation of this algorithm is demonstrated, and emphasizes the practical benefits of this new method over the previous state of the art.

We consider the problem of sparse interpolation of an approximate multivariate black-box polynomial in floating-point arithmetic. That is, both the inputs and outputs of the black-box polynomial have some error, and all numbers are represented in standard, fixed-precision, floating point arithmetic. By interpolating the black box evaluated at random primitive roots of unity, we give efficient and numerically robust solutions. We note the similarity between the exact Ben-Or/Tiwari sparse interpolation algorithm and the classical Prony's method for interpolating a sum of exponential functions, and exploit the generalized eigenvalue reformulation of Prony's method. We analyze the numerical stability of our algorithms and the sensitivity of the solutions, as well as the expected conditioning achieved through randomization. Finally, we demonstrate the effectiveness of our techniques in practice through numerical experiments and applications.

New algorithms are presented for computing the dispersion set of two polynomials over Q and for shiftless factorization. Together with a summability criterion by Abramov, these are applied to get a polynomial-time algorithm for indefinite rational summation, using a sparse representation of the output.

As a sub-procedure our algorithm executes the Berlekamp/Massey algorithm on a sequence of large integers or polynomials. We give a fraction-free version of the Berlekamp/Massey algorithm, which does not require rational numbers or functions and GCD operations on the arising numerators and denominators. The relationship between the solution of Toeplitz systems, Padé approximations, and the Euclidean algorithm is classical. Fraction-free versions [3] can be obtained from the subresultant PRS algorithm [2]. Dornstetter [6] gives an interpretation of the Berlekamp/Massey algorithm as a partial extended Euclidean algorithm. We map the subresultant PRS algorithm onto Dornstetter's formulation. We note that the Berlekamp/Massey algorithm is more efficient than the classical extended Euclidean algorithm.

A new algorithm is presented for factoring bivariate approximate polynomials over C[x, y]. Given a particular polynomial, the method constructs a nearby composite polynomial, if one exists, and its irreducible factors. Subject to a conjecture, the time to produce the factors is polynomial in the degree of the problem. This method has been implemented in Maple, and has been demonstrated to be efficient and numerically robust.

A probabilistic algorithm is presented to find the determinant of a nonsingular, integer matrix. For a matrix A/spl isin/Z/sup n/spl times/n/ the algorithm requires O(n/sup 3.5/(log n)/sup 4.5/) bit operations (assuming for now that entries in A have constant size) using standard matrix and integer arithmetic. Using asymptotically fast matrix arithmetic, a variant is described which requires O(n/sup 2+/spl theta//2//spl middot/log/sup 2/nloglogn) bit operations, where n/spl times/n matrices can be multiplied with O(n/sup /spl theta//) operations. The determinant is found by computing the Smith form of the integer matrix an extremely useful canonical form in itself. Our algorithm is probabilistic of the Monte Carlo type. That is, it assumes a source of random bits and on any invocation of the algorithm there is a small probability of error.

The (4 0,) (1:) or fadoretl fornl.For exau~ple~ a. demise polynomial of degree IL woultl take approsimately 2r1 operat.ionst.o cva1uat.e in eitlwr espandcd or factin form.-4 presentatiori itS two cornlx~sit.ionfac:tor.r;,however: ~voldtl t,akr l>etw:rn 4&i and II aritlm&ic operations.main results of this pilpcr illC ail iterative nictliod t,o conlput~e il.decomposit.ion of a giveu itpprosillla~.e pOl~IlOIlliill, giveIl a starting point.

Randomized black b o x algorithms provide a very ecient means for solving sparse linear systems over arbitrary elds.However, when these probabilistic algorithms fail, it is not revealed whether no solution exists or whether the algorithm simply made unlucky random choices.Here we give a n ecient algorithm to compute a certi cate of inconsistency for a black b o x linear system over a eld.Our method requires a black box for the transpose of the matrix.The cost of producing the certi cate is shown to be about the same as that of solving the system in the black b o x model, while the cost of applying a given certi cate to prove inconsistency is much smaller.We also give an ecient algorithm for certifying that a sparse Diophantine linear system of integer equations has no integer solutions, even when it may h a v e rational solutions.

We present new, efficient algorithms for some fundamental computations with finite-dimensional (but not necessarily commutative)associative algebras.For a semisimple associative algebra ~over a finite field or number field F, we show how to compute a basis for the centre of X and the complete Wedderburn decomposition of 2! as a direct sum of simple algebras.If ~is given by a generating set of matrices in F~xm then our algorithm requires about 0(rn3) operations in F, plus the cost of factoring a polynomial in F[z] of degree O(m), and the cost of generating a small number of random elements from fl.We also show how to compute a complete set of orthogonal primitive idempotents in any associative algebra over a finite field.

We present two new probabilistic algorithms for computing the Smith normal form of an A 2 Z m\\Thetan . The first requires an expected number of O(m 2 n \\Delta M(m log kAk)) bit operations (ignoring logarithmic factors) and is of the Las Vegas type; that is, it never produces an incorrect answer. Here kAk = max ij jA ij j and M(l) bit operations are sufficient to multiply two l-bit integers (M(l) = l 2 using standard arithmetic) . This improves on the previously best known (deterministic) algorithm of Hafner and McCurley, which requires about O(m 3 n log kAk \\Delta M(m log kAk)) bit operations. We also present an even faster, more space efficient algorithm which requires an expected number of O((m 3 n log kAk + m 3 log 2 kAk) \\Delta log(1=ffl)) bit operations using standard integer arithmetic. This algorithm is of the Monte Carlo type: it returns the correct result with probability at least 1 \\Gamma ffl for a user specified tolerance ffl ? 0. This algorithm also require...

A Las-Vegas-type probabilistic algorithm is presented for finding the Frobenius canonical form of an $n \times n$ matrix T over any field $\mathfrak{K}$. The algorithm requires $O^{\sim}(\operatorname{MM}(n)) = \operatorname{MM}(n) \cdot (\log n)^{O(1)}$ operations in $\mathfrak{K}$, where $O(\operatorname{MM}(n))$ operations in K are sufficient to multiply two $n \times n$ matrices over $\mathfrak{K}$. This nearly matches the lower bound of $\Omega (\operatorname{MM}(n))$ operations in $\mathfrak{K}$ for this problem and improves on the $O(n^{4})$ operations in $\mathfrak{K}$ required by the previous best-known algorithms. A fast parallel implementation of the algorithm is also demonstrated for the Frobenius form, which is processor-efficient on a PRAM. As an application we give an algorithm to evaluate a polynomial $g \in {\mathfrak{K}}[x]$ at T which requires only $O^{\sim}(\operatorname{MM}(n))$ operations in $\mathfrak{K}$ when $\deg g \leq n^{2}$. Other applications include sequential and parallel algorithms for computing the minimal and characteristic polynomials of a matrix, the rational Jordan form of a matrix (for testing whether two matrices are similar), and matrix powering, which are substantially faster than those previously known.

A Monte Carlo type probabilistic algorithm is presented for finding the Frobenius rational form F 2 Z n\\Thetan of any A 2 Z n\\Thetan which requires an expected number of O(n 4 (log n + kAk) 2 ) bit operations using standard integer and matrix arithmetic (where kAk is the largest absolute value of any entry of A). This improves dramatically on the fastest previously known algorithm, which requires O(n 6 log kAk) bit operations using fast integer arithmetic. We also give a Las Vegas type probabilistic algorithm which finds the Frobenius form F and a transition matrix U 2 Q n\\Thetan such that U \\Gamma1 AU = F and requires an expected number of O(n 5 (log n + log kAk) 5=2 bit operations. Finally, a Las Vegas algorithm for computing the rational Jordan form of an integer matrix is shown, which requires about the same number of bit operations as our algorithm to find the Frobenius form, plus the time required to factor the characteristic polynomial of that matrix. 1. I...

A Las Vegas type probabilistic algorithm is presented for computing the Frobenius normal form of an n*n matrix T over any field K. The algorithm requires O/sup approximately /(MM(n))=MM(n)/sup ./(log n)/sup O(1)/ operations in K, where O(MM(n)) operations in K are sufficient to multiply two n*n matrices over K. This nearly matches the lower bound of Omega (MM(n)) operations in K for this problem, and improves on the O(n/sup 4/) operations in K required by the previously best known algorithm. The author applies the algorithm to evaluate a polynominal g in K(x) at T with /sup approximately /(MM(n)) operations in K when deg g>