# Computing a polynomial product over roots of unity

I'm trying to compute the coefficients of the following polynomial, where $\omega$ is a primitive $p$-th root of unity, for $p$ prime:

$$a(x) = \prod_{i=0}^{p-1} f(\omega^ix).$$

It turns out that the $i$-th coefficient is always an integer, and non-zero only when $i$ is a multiple of $p$. So it seems to me like there should be an elementary expression for $a$.

So far I've got this expression for the $i$-th coefficient: $$a_i = x^i\sum_{k_0 + \ldots + k_{p-1} = i}f_{k_0} \cdots f_{k_{p-1}} \omega^{k_1 + 2k_2 + \ldots + (p-1)k_{p-1}}$$

where each $k_i$ is non-negative and bounded by the degree of $f$.

Clearly the roots of unity all cancel out somehow, but I can't figure out how to get a 'nice' expression out. Any suggestions?

-
Do you know the roots of $f$? If $f(x)=c\prod_{i < d}(x-\alpha_i)$, then $a$ is easily seen to be $(-1)^{d(p+1)}c^p\prod_{i < d}(x^p-\alpha_i^p)$. –  Emil Jeřábek Mar 28 '11 at 16:47
No, unfortunately I don't know the roots of $f$. Thanks anyway though, that's interesting to know. –  Peter Scholl Mar 28 '11 at 17:18
The last (i.e. degree zero) coefficient is just the $p$-th power of the last coefficient of $f$. Other than that, I suspect that there's no "nice" expression for the coefficients of $a$: for example, already for the first nontrivial coefficient (the one in degree $p\cdot\deg f$), the expression involves all of the coefficients of $f$. Computing this coefficient is related to Newton's identities (see en.wikipedia.org/wiki/Newton%27s_identities ). –  Marco Golla Mar 28 '11 at 17:23
We may assume that $f(0)=1$. If $g(x)=\log f(x) = g_1x+g_2x^2+\cdots$ then $\log a(x) = p(g_px^p+g_{2p}x^{2p}+\cdots)$. We can then compute $a(x)$ by exponentiation. This will give a formula for $a(x)$ not involving roots of unity, but it won't be pretty. On the other hand, there are fast algorithms for computing log and exp of a power series. –  Richard Stanley Mar 28 '11 at 17:39
"The $i$th coefficient is always an integer, and non-zero only when $i$ is a multiple of $p$". I guess $f\in\mathbf{Z}[x]$? If so then the coefficients are rational, by Galois theory, and integers because they're algebraic integers, and vanish if $i$ isn't a multiple of $p$ because if $F$ is the product then $F(\omega x)=F(x)$. However, "seems to me like there should be an elementary expression"---doesn't really follow as far as I can see. Take a random recurrence relation like $a(n+2)=n^2*a(n+1)+a(n)^2$ with $a(1)=a(2)=38$. Then all $a(n)$ are clearly integers but why expect a nice formula? –  Kevin Buzzard Mar 28 '11 at 18:43

Notice that $$y^p-x^p=\prod_{i=0}^{p-1} (y-\omega^i x)\ .$$ Therefore $$\prod_{i=0}^{p-1} f(x\omega^i)= Res_y (y^p-x^p, f(y))\ .$$ Where $Res_y$ is the resultant of the polynomials in $y$. The above is just a particular case of the so called Poisson product formula. You can then compute the resultant using for instance Sylvester's determinant formula. See the review http://mate.dm.uba.ar/~alidick/papers/chapter1cd.pdf by Cattani and Dickenstein for a nice introduction to resultants and their properties.

-
That's a nice result, but I think computing the resultant here is quite slow when $f$ is large, which I need it to be. What I'm really looking for (perhaps I should have made it more clear) is an efficient way to recover $a$, rather than just an expression. –  Peter Scholl Mar 29 '11 at 12:36
@Peter: It is unclear to me what you really want. Is there (1) a specific f for which you would like to compute a(x), or do you want (2) a formula for a in terms of a generic f? If it is (1) then the resultant calculation is that of a determinant of order p+deg(f) with numerical coefficients and a lot of zeros, shouldn't be too bad with a computer algebra software. If it is (2), I don't think you will find a formula for a in terms of undetermined coefficients of f which is simpler that the expansion of the determinant giving the resultant. –  Abdelmalek Abdesselam Mar 29 '11 at 15:10
@Peter : there is a division-like algorithm to compute the resultant (well, the division algorithm only works for polynomials over a field, but there is a notion of subresultant polynomial, see mathworld.wolfram.com/Subresultant.html which circumvents this difficulty). In general I doubt that there will be a simpler formula than the one given by Abdelmalek, unless you have some more information on f. –  François Brunault Mar 30 '11 at 7:49

I don't think that $p$ being prime makes any difference. The later thoughts below suggest a Lagrange interpolation method which is perhaps the same as the resultant method mentioned by Abdelmalek Abdesselam.

Let $f(x)=\sum_{j=0}^nc_jx^j.$ One might require $c_n=1$ or $c_0=1$ but it is perhaps nicer not to. Then setting $u=x^n,$ $\prod_{i=0}^{n-1} f(\omega^ix)=F(u)=\sum_{j=0}^nC_j u^j$. One can say that

• $C_j$ is a polynomial of degree $n$ in the coefficients $c_0,\cdots,c_n$ where each term has total degree $n$
• $C_j$ has a term $\pm(c_j)^n$ and no term $c_j^{n-1}$.
• $C_{n-j}$ is $C_j$ with $c_k$ replaced by $c_{n-k}$
• the roots of $F$ are the $n$th powers of the roots of $f$.

Here $c_n$ is a constant and the other $c_i$ are symmetric polynomials of the $n$ roots $\alpha_i$ of $f$. The $\alpha_i$ can be thought of as formal variables. Then $c_0,\cdots,c_{n-1}$ are also a basis for the ring of all symmetric polynomials in those variables ($\frac{1}{c_n}$ times a usual basis) . There are other bases for this ring such as $\sigma_k=\sum_{i=1}^n\alpha_i^k$ and the sum of all the terms $\alpha_1^{m_1}\cdots\alpha_n^{m_n}$ with $m_1+\cdots+m_n=k$. Transforming between these bases (more generally, expressing a given symmetric polynomial in terms of them) is a major topic of invariant theory.

In this case, we want to express the $C_i$, which are certain symmetric polynomials in the $\alpha_i^n,$ in terms of the values $c_i$. This must be a well known case. At any rate here are some results:

For $n=4,$

$$C_0=c_0^4$$

$$C_1=(4c_0^3c_4-2c_0^2c_2^2)-(4c_0^2c_3)c_1+(4c_0c_2)c_1^2-c_1^4$$

$$C_2=(6c_0^2c_4^2-8c_0c_1c_3c_4+2c_1^2c_3^2)+(4c_0c_3^2+4c_1^2c_4)c_2-(4c_0c_4+4c_1c_3)c_2^2+c_2^4$$

$$C_3=(4c_4^3c_0-2c_4^2c_2^2)-(4c_4^2c_1)c_3+4(c_4c_2)c_3^2-c_3^4$$

$$C_4=c_4^4$$

while for $n=5$ we have the following (with the others obtainable by symmetry)

$$C_0=c_0^5$$

$$C_1=5c_0^4c_5-(5c_0^3c_4-5c_0^2c_2^2)c_1+(5c_0^2c_3)c_1^2-5c_0^3c_2c_3-(5c_0c_2)c_1^3+c_1^5$$

$${\small C_2=(10c_0^3c_5^2-15c_0^2c_1c_4c_5+5c_0c_1^2c_4^2+5c_0^2c_3^2c_4+10c_0c_1^2c_3c_5-5c_0c_1c_3^3-5c_1^3c_3c_4)}$$ $${ \small+(5c_0^2c_4^2-15c_0^2c_3c_5+5c_1^2c_3^2-5c_0c_1c_3c_4-5c_1^3c_5)c_2}$$ $${\small +(5c_1^2c_4+10c_0c_1c_5+5c_0c_3^2)c_2^2-(5c_0c_4+5c_1c_3)c_2^3+c_2^5}$$

later thoughts In general one could consider the problem of producing from a polynomial $f(x)=\sum_{j=0}^nc_jx^j$ a polynomial $F(u)=\sum_{j=0}^nC_j u^j$ whose roots are the $m$th powers $\alpha_i^m$ of the $n$ (unknown) roots of $f$. The solution is the polynomial $\prod_{i=0}^{m-1} f(\tau^ix)$ where now $\tau$ is a primitive $m$th root of unity. This process might spread out the roots. In the case $m=2$ one has (with $q$ repeated application) $\alpha_i^{2^q}$ and the Dandelin–Graeffe method for finding the roots of a univariate polynomial. Splitting $f$ into its even and odd parts speeds up the computation. The method was also discovered by Nikolai Ivanovich Lobachevsky and the linked article suggests that his book Algebra ili Ichislenie Konechnykh Velichin discussed the general product $\prod_{i=0}^{m-1} f(\omega^ix)$. Perhaps the appropriate manipulations of symmetric polynomials are discussed there.

Since a polynomial of degree $n$ is determined by its values at $n$ points (plus its leading coeffcient) one has the following method (which I doubt is new): Let $\zeta$ be a primitive $mn$th root of unity and $\omega$ a primitive $m$th root of unity. Then the desired polynomial $F$ satisfies $F(\zeta^j)=f(\omega^j)$ for $0 \le j \le n-1$. Now Lagrange Interpolation can be used. In this case it might be particularly simple.

-