Roots multiplicity without companion matrices

TL;DR

This paper introduces a new polynomial interpolation method for roots' multiplicities that avoids companion matrices, improving the efficiency of the square-free decomposition algorithm in both theory and practice.

Contribution

It presents a novel approach to roots' multiplicity interpolation that enhances the Guersenzvaig--Szechtman algorithm without relying on companion matrices.

Findings

More efficient square-free decomposition algorithm

Effective polynomial interpolation method for roots' multiplicities

Practical improvements demonstrated in computational experiments

Abstract

We show a method for constructing a polynomial interpolating roots' multiplicities of another polynomial, that does not use companion matrices. This leads to a modification to Guersenzvaig--Szechtman square-free decomposition algorithm that is more efficient both in theory and in practice.

Tables

Table 1. Table 1. Running times (in seconds) of Magma implementations of Eqs. ( A ) and ( B ). For each degree, both formulas were evaluated on the same set of 10 randomly generated polynomials.

degree	Eq. (A)		Eq. (B)
10	0.	010	0.	000
20	0.	020	0.	010
50	0.	170	0.	140
100	1.	140	0.	600
200	10.	690	3.	170
500	330.	900	32.	230

Equations

$$f=P_1\cdot P_2^2\cdots P_m^m,$$

$$f=\prod _{\alpha \in S(f)}(X-\alpha {)}^{m(\alpha )}.$$

$$P_k=\gcd (M_f-k,r)\text{for }k\in \{1,\dots ,m\}.$$

$$r=\prod _{\alpha \in S(f)}(X-\alpha )=\frac{f}{\gcd (f,f^{\prime })}.$$

$$C_r:=\left(\begin{array}{cccc}0& & & -r_0\\ 1& \ddots & & ⋮\\ & \ddots & 0& -r_{s-2}\\ & & 1& -r_{s-1}\end{array}\right)$$

$$P:=\frac{f^{\prime }}{\gcd (f,f^{\prime })}.$$

$$r^{\prime }\cdot g+r\cdot h=1\text{and}\mathrm{deg}g<\mathrm{deg}r,\mathrm{deg}h<\mathrm{deg}{r}^{\prime }.$$

$$[M_f]=P(C_r)\cdot [g].$$

$$P(\alpha )=m(\alpha )\cdot r^{\prime }(\alpha ).$$

$$r^{\prime }=\sum _{\alpha \in S(f)}\prod _{\begin{array}{c}\beta \in S(f)\\ \beta \ne \alpha \end{array}}(X-\beta ).$$

$$r^{\prime }(\alpha )=\prod _{\begin{array}{c}\beta \in S(f)\\ \beta \ne \alpha \end{array}}(\alpha -\beta ).$$

$$f^{\prime}=\sum_{\alpha\in S(f)}\Bigl{(}m(\alpha)\cdot(X-\alpha)^{m(\alpha)-1}\cdot\prod_{\begin{subarray}{c}\beta\in S(f)\\ \beta\neq\alpha\end{subarray}}(X-\beta)^{m(\beta)}\Bigr{)}.$$

$$P$$

$$\displaystyle=\frac{\sum\limits_{\alpha\in S(f)}\Bigl{(}m(\alpha)\cdot(X-\alpha)^{m(\alpha)-1}\cdot\prod\limits_{\beta\neq\alpha}(X-\beta)^{m(\beta)}\Bigr{)}}{\prod\limits_{\alpha\in S(f)}(X-\alpha)^{m(\alpha)-1}}$$

$$\displaystyle=\sum_{\alpha\in S(f)}\Bigl{(}m(\alpha)\cdot\prod_{\begin{subarray}{c}\beta\in S(f)\\ \beta\neq\alpha\end{subarray}}(X-\beta)\Bigr{)}.$$

$$P(\alpha )=m(\alpha )\cdot \prod _{\begin{array}{c}\beta \in S(f)\\ \beta \ne \alpha \end{array}}(\alpha -\beta )$$

$$M_{f}=\bigl{(}P\cdot g\bmod r\bigr{)}.$$

$$\mu (\alpha )=P(\alpha )\cdot g(\alpha )=m(\alpha )\cdot r^{\prime }(\alpha )\cdot g(\alpha )=m(\alpha )\cdot (1-r\cdot h)(\alpha )=m(\alpha ).$$

$$P:=\frac{f^{\prime }}{\gcd (f,f^{\prime })},r:=\frac{f}{\gcd (f,f^{\prime })}$$

$$P_{k}:=\gcd\bigl{(}M_{f}-k,r\bigr{)}.$$

Full text

The problem of computing the square-free decomposition of a polynomial is well established in the realm of computational algebra. There are efficient algorithms for this task developed nearly half a century ago by: R. Tobey [6], E. Horowitz [3], D. Musser [5] and D. Yun [7]. Recently, N. Guersenzvaig and F. Szechtman invented a completely new algorithm (see [2]). They associate to a given polynomial $f$ its roots-multiplicity polynomial $M_{f}$. The authors presented a formula, for constructing $M_{f}$, based on the companion matrix of the radical of $f$. The aim of this note is to show that $M_{f}$ can be obtained much more efficiently without using the companion matrix.

For readers convenience (and to some discomfort of the author), in this paper we use the same notation as N. Guersenzvaig and F. Szechtman did. Let $F$ be a fixed field of characteristic [math]. Given a monic polynomial $f\in F[X]$, the square-free decomposition of $f$ is an expression

[TABLE]

where each $P_{k}\in F[X]$ for $k\in\{1,\dotsc,m\}$ is monic and square-free. Let $S(f)$ be the set of all the roots of $f$ (in some algebraically closed field). For every root $\alpha\in S(f)$, denote by $m(\alpha)$ the multiplicity of $\alpha$, so that

[TABLE]

By Lagrange interpolation formula, there is a unique polynomial $M_{f}$ of minimal degree such that $M_{f}(\alpha)=m(\alpha)$ for every root $\alpha$ of $f$. Once the polynomial $M_{f}$ is found, the square free factors $P_{1},\dotsc,P_{m}$ of $f$ can be computed by the formula [2, Eq. (1.5)]:

[TABLE]

Thus, all we need is an efficient method for constructing $M_{f}$.

Guersenzvaig and Szechtman proposed the following procedure. The radical $r:=\operatorname{rad}f$ is a polynomial

[TABLE]

The radical is square-free, hence its degree, denote it by $s:=\deg r$, equals the cardinality of $S(f)=S(r)$. Suppose that $r$ has a form $r=r_{0}+r_{1}X+\dotsb+r_{s-1}X^{s-1}+X^{s}$. Let

[TABLE]

be a companion matrix of $r$ and set

[TABLE]

The fact that the radical is square-free implies that it is relatively prime to its derivative. Therefore, by Bézout identity, there are polynomials $g,h\in F[X]$ such that

[TABLE]

For any polynomial $p$, $\deg p<s$, by $[p]$ denote a column vector of its coefficients, zero-appended to length $s$, if needed.

[TABLE]

We claim that $M_{f}$ can be constructed more efficiently. To this end we need:

[TABLE]

[TABLE]

It follows that for every root $\alpha$ of $f$ we have

[TABLE]

Next, compute the derivative of $f$:

[TABLE]

It follows that the polynomial $P$ can be expressed in a form

[TABLE]

Consequently, evaluating $P$ at a root $\alpha$, we have

[TABLE]

and the thesis follows from Eq. (3). ∎

[TABLE]

[TABLE]

This shows that these two polynomials are identical. ∎

Let us now recall Guersenzvaig–Szechtman algorithm.

1. (a)

Construct polynomials

[TABLE]

and find polynomials $g,h$ satisfying Eq. (2); 2. (b)

build the polynomial $M_{f}$ using either formula (A) or formula (B); 3. (c)

for $k=1,2,\dotsc$ and as long as $\sum_{j\leq k}j\cdot\deg P_{j}<\deg f$, compute $k$-th square-free factor of $f$:

[TABLE]

Complexity analysis

As it is a common practice, let $M(n)$ denote the time complexity of computing a product of two polynomials of degree $n$. For long polynomial multiplication $M(n)=O(n^{2})$, while for fft-based multiplication $M(n)=O(n\cdot\lg n)$. (Here $\lg=\log_{2}$ stands for logarithm of base $2$.) Recall that $s=\deg r$ and let $n:=\deg f$. It is clear that $s=O(n)$. If one uses fast gcd, then step (a) of the algorithm has time complexity O\bigl{(}M(n)\lg n\bigr{)}. The same applies to every iteration in step (c). The number of iterations equals $m=O(n^{0.5})$, the number of square-free factors of $f$. Hence the complexity of step (c) is O\bigl{(}n^{0.5}\cdot M(n)\cdot\lg n\bigr{)}=O\bigl{(}n^{1.5}\cdot(\lg n)^{2}\bigr{)}, when using fast polynomial multiplication and $O(n^{2.5}\lg n)$ for long multiplication. It remains to analyze the time complexity of step (b), the one for which we propose a new formula.

Observe that in order to construct $M_{f}$ with formula (A), one needs to evaluate the polynomial $P$ on the matrix $C_{r}$. Using Horner scheme, this requires $s$ products of $s\times s$ square matrices. Thus, with standard matrix arithmetic, Eq. (A) needs $s+s^{4}$ scalar multiplications. Using fast matrix multiplication (see e.g. [4, 8]), the asymptotic complexity (the actual number of scalar products is hard to count) of that formula is about $O(s^{3.37})=O(n^{3.37})$. Therefore, if formula (A) is in use, step (b) asymptotically dominates the running time of the whole algorithm.

On the contrary, formula (B) needs only $s+M(s)$ scalar products and so it has the asymptotic complexity of $O(s\lg s)=O(n\lg n)$, if one uses fast polynomial multiplication (respectively $O(s^{2})=O(n^{2})$ for long multiplication). Consequently, with formula (B), step (b) no longer dominates the time complexity of the algorithm.

In order to compare both formulas in real-life situations, they were implemented in computer algebra system Magma [1] and evaluated on random polynomials of different degrees. Table 1 presents the running times. The code can be downloaded from author’s website http://z2.math.us.edu.pl/perry/papersen.html

Conclusion

The modification of Guersenzvaig–Szechtman algorithm presented in this note leads to a procedure that is faster than the original one both theoretically and empirically.