Arithmetic on algebraic numbers

I'd like to perform some elementary operations on algebraic numbers.

p1 = x^5 + 6*x^4 - 42*x^3 - 142*x^2 + 467*x + 422;
p2 = Expand[p1 /. {x -> ((x - 1)^2)}];
r1 = Root[{p2 /. x -> # &, 1 + 2.14850301089246970680618135585609253400250*I}];
r2 = Root[{p2 /. x -> # &, 1 + 2.87513448574410025133873802142153923543448*I}];
Timing[MinimalPolynomial[r1-r2]]

I'm trying to compare Mathematica's behavior to what I'm currently working on with Sage. In particular I'd like to see whether division r1/r2 takes much longer than r1+r2, r1-r2, r1*r2.

But it turns out that Mathematica apparently isn't able to perform elementary arithmetic operations on these algebraic numbers at all. The above example fails with

MinimalPolynomial::nalg: "-Root[…]+Root[…] is not an explicit algebraic number."

Is there some formulation I can use to perform such computations? Or is Mathematica simply not up to this task?

• Welcome to Mathematica.SE! I suggest that: 1) You take the introductory Tour now! 2) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! 3) As you receive help, try to give it too, by answering questions in your area of expertise. – bbgodfrey Mar 4 '15 at 13:42

Root[{f, approx}] format was meant for representing roots of non-polynomial functions. Mathematica will never produce algebraic numbers (of degree <= \$MaxRootDegree) represented in this way. If you manually enter Root[{f, approx}] with a polynomial f and a correct approximation the result is a valid exact number, but it is not recognized as an algebraic number. I will make RootReduce check for such objects and convert them to the standard algebraic number representation -- Root[poly, root number]. For now you can automate the method suggested by Daniel.

In:= rootFromApprox[poly_, x_, approx_]:=
Select[RootReduce[x/.Solve[poly==0, x, Cubics->False, Quartics->False]],
#-approx==0&][]

In:= p1=x^5+6*x^4-42*x^3-142*x^2+467*x+422;
p2=Expand[p1/.{x->((x-1)^2)}];
r1=rootFromApprox[p2, x, 1+2.14850301089246970680618135585609253400250*I];
r2=rootFromApprox[p2, x, 1+2.87513448574410025133873802142153923543448*I];
Timing[Short[MinimalPolynomial[#]]]&/@{r1+r2, r1-r2, r1 r2, r1/r2}

Out= {{0.124801,-3568379502241854768+<<58>>+#1^40&},
{0.140401,<<34>>+48 #1^38+#1^40&},
{0.140401,66045000696445844586496+82370731205679648866304 #1+<<58>>+#1^40&},
{0.327602,<<118>>+257988283970491580416 #1^80&}}

I do not know what is going awry and have sent a question to the relevant person here. I can suggest an alternative route. Form the Root objects using say Solve, figure out which are the ones you want, and use those.

p1 = x^5 + 6*x^4 - 42*x^3 - 142*x^2 + 467*x + 422;
p2 = Expand[p1 /. {x -> ((x - 1)^2)}];
rts = RootReduce[x /. Solve[p2 == 0, x]];
N[rts]

(* Out= {1. - 2.87513448574 I, 1. + 2.87513448574 I,
1. - 2.14850301089 I, 1. + 2.14850301089 I, 1. - 0.876768035257 I,
1. + 0.876768035257 I, -0.823985071296, 2.8239850713, \
-1.07948651024, 3.07948651024} *)

So it is the fourth and second that we want.

{rt1, rt2} = rts[[{4, 2}]]
Timing[MinimalPolynomial[rt1 - rt2]]

(* Out= {Root[
712 - 172 #1 - 802 #1^2 + 952 #1^3 - 142 #1^4 - 336 #1^5 +
336 #1^6 - 168 #1^7 + 51 #1^8 - 10 #1^9 + #1^10 &, 8],
Root[712 - 172 #1 - 802 #1^2 + 952 #1^3 - 142 #1^4 - 336 #1^5 +
336 #1^6 - 168 #1^7 + 51 #1^8 - 10 #1^9 + #1^10 &, 6]}

Out= {0.050724,
515591535280592 - 6604303796727744 #1^2 - 19213629465035648 #1^4 -
9238448577096832 #1^6 - 246985961012832 #1^8 -
408212448315616 #1^10 - 179542482936184 #1^12 -
20639353188832 #1^14 - 6030161483696 #1^16 - 874984589616 #1^18 -
187565974272 #1^20 - 33966863800 #1^22 - 4586283639 #1^24 -
530873932 #1^26 - 40291328 #1^28 - 2037008 #1^30 + 5146 #1^32 +
10684 #1^34 + 924 #1^36 + 48 #1^38 + #1^40 &} *)