Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

How does Mathematica's RootOfUnityQ function work?

That is, how does Mathematica know if a number is a root of unity?

Example:

Let $x = \frac{1-i \sqrt{2+\sqrt{5}}}{1+i \sqrt{2+\sqrt{5}}}$.

Then

In[1]:= RootOfUnityQ[x]
Out[1]= False

Wrong Answers:

This bahavior is not the same as

In[2]:= Reduce[Abs[x] == 1 && Simplify[Element[Arg[x]/(2 Pi), Rationals]]]
Out[2]= ArcTan[(2 Sqrt[2 + Sqrt[5]])/(1 + Sqrt[5])]/\[Pi] \[Element] Rationals
share|improve this question
    
This doesn't answer the question of how Mathematica does it, but we discussed some potential methods at math.stackexchange.com/questions/68159 –  David Speyer Aug 14 '12 at 23:59
add comment

2 Answers 2

As long as we're giving methods that mimic the behavior, here is a quick implementation based on Kronecker's theorem.

RoUQ[u_] :=
 If[! (Abs[N[u]] == 1), False,
  If[! AlgebraicIntegerQ[u], False,
   (f = MinimalPolynomial[u, x];
    n = Exponent[f, x];
    cf = CoefficientList[f, x]/Coefficient[f, x^n]; 
    M = Table[
      If[j == i + 1, 1,
       If[i == n, -cf[[j]], 0]], {i, 1, n}, {j, 1, n}];
    lambda = First[Eigenvalues[M, 1]];
    Abs[N[lambda]] == 1)]] 

The sole reason for building the matrix $M$ is that Eigenvalues[M,1] is guaranteed to return the largest eigenvalue of $M$ (in absolute value), while I don't know how to tell FindRoot that it must give me the largest root.


J.M. below suggests just finding all the roots and taking the largest one. Based on a tiny amount of testing, I think that this is slower. Here is what I did:

(* Both methods will be tested on the same set of polynomials *)
data=Table[Sum[RandomInteger[100] x^j, {j,0,10}], {i,1,10^4}];

(* Find all roots and extract the largest *)
BigRoot1[f_]:=Max[Abs[x] /. NSolve[f==0,x]]

(* Find the largest eigenvalue of the companion matrix. *)
BigRoot2[f_]:=
(n = Exponent[f, x];
    cf = CoefficientList[f, x]/Coefficient[f, x^n]; 
    M = Table[
      If[j == i + 1, 1,
       If[i == n, -cf[[j]], 0.0]], {i, 1, n}, {j, 1, n}];
    lambda = First[Eigenvalues[M, 1]];
    N[Abs[lambda]])                                                               

Timing[Map[BigRoot1, data]][[1]]

10.9893

Timing[Map[BigRoot2, data]][[1]]

3.00254

I wouldn't take these results too seriously, because I'm sure both implementations do the list processing inefficiently, but it suggests that I'm not crazy to use Eigenvalues[].

share|improve this answer
    
Well, you could use x /. First @ Solve[p[x] == 0, x] to generate the roots, and then Pick[roots, #, Max[#]] &[Abs[roots]] to get the largest of them all in magnitude. But using the companion matrix is cool anyway. :) –  J. M. Aug 15 '12 at 23:37
    
@J.M. It strikes me as wasteful to compute all those roots which will be discarded. I've added some tests above to support this claim. –  David Speyer Aug 16 '12 at 0:14
    
Hmm, right. Then, I suppose if one wants to use FindRoot[], there's always the option of trying to find the smallest root of $x^n p\left(\frac1{x}\right)$ and reciprocating that... the only reason I was somewhat lukewarm to using the companion matrix is that it looks uncomfortable to use $O(n^2)$ storage and $O(n^3)$ effort for just a single root... –  J. M. Aug 16 '12 at 1:45
    
@J.M. Yeah, memory usage is a good point. I guess the next question then is how often FindRoot[], starting at 0, fails to find the smallest root. A good first test would be p[x]=x^10-x^9+x^7-x^6+x^5-x^4+x^3-x+1, or perhaps this after the change of variable x /. E^(2 Pi I RandomReal[]) x. This has two roots at 0.85 and 1.176 and the other 8 are all on the unit circle with roughly equal spacing, so a random walk from 0 is very likely to hit one of the roots on the circle before the small root. I imagine FindRoot[] does something smarter than a random walk, though. –  David Speyer Aug 16 '12 at 12:13
add comment

If Abs[x]==1 and Element[Arg[x]/(2 Pi), Rationals], then the number is a root of unity

$$ \arg[x]= 2 \pi \frac{p}{q}$$ $p,q\in \mathbb{Z}$

and $$ |x|=1 $$

$$ \to x^q=1 $$

because
$ |x^q|=|x|^q=1^q=1 $ and
$ \arg[x^q] = \arg[x] q = 2 \pi p $

So perhaps this mimics the behaviour

ru = TrueQ[Abs[#] == 1 && Simplify@Element[Arg[#]/(2 Pi), Rationals] ]&;
share|improve this answer
    
It does not mimic the behavior. I updated the question to discuss your answer. –  Tyson Williams Aug 1 '12 at 19:48
    
@TysonWilliams in your counterexample, it's just a matter of adding TrueQ. There are probably other counterexamples where this doesn't work, but please post them. I'm editing my answer –  Rojo Aug 1 '12 at 20:08
    
@R.M I am not looking for (numerically) approximate answers. I want to know with certainty. –  Tyson Williams Aug 1 '12 at 20:10
    
Ok, suppose that this is the definition of RootOfUnityQ. The purpose of my question was to be able to verify Mathematica's answer for myself. So my question now becoems this: how does Mathematica compute ArcTan[Sqrt[1/2 (1+Sqrt[5])]]/\[Pi]\[Element]Rationals? –  Tyson Williams Aug 1 '12 at 21:34
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.