Solving Cyclotomic Polynomials by Radical Expressions
Andreas Weber, Michael Keckeisen, Essam Abdel-Rahman
E-mail: weber@cs.uni-bonn.de
Introduction
Niels Henrik Abel and Evariste Galois showed that polynomial equations of degree higher than four cannot be solved by radical expressions in general. As Galois stated in his work, radical solutions exist if and only if the Galois group of the polynomial is solvable.
Since the the Galois group of a cyclotomic polynomial is Abelian, its Galois group is solvable, and so its solutions can be expressed by radical expressions.
Example: the 5th cyclotomic polynomial and it's solutions:
| > |
cycl:=numtheory[cyclotomic](5,x); |
| > |
for i from 1 to 4 do t[i]:=[Re(s[i]),Im(s[i])]; u[i]:=[Re(s[i]),Im(s[i]), ZETA^i] od: |
| > |
p1:=plot([[1,0],t[1],t[2],t[3],t[4],[1,0]],style=line,title=`The 4 Solutions of the 5th cyclotomic polynomial:`): |
| > |
p2:=plots[implicitplot](x^2+y^2=1,x=-1..1,y=-1..1,scaling=constrained): |
| > |
p3:=plots[textplot]([[1,0,``],u[1],u[2],u[3],u[4],[1,0,``]]): |
| > |
plots[display]([p1,p2,p3]); |
In the last years algorithmic methods to compute solutions by radicals for the cases in which the Galois group of the polynomial is solvable have obtained renewed attention [HM01]. However, for the special case of cyclotomic polynomials a very good algorithm was already published in the year 1801, when the masterwork Disquisitiones Arithmeticae ([Gau86]) of Carl Friedrich Gauss appeared. In this book Gauss described an algorithm to compute radical expressions for primitive roots of unity. This work in some sense is a generalization of the famous result of Gauss already obtained in 1797 showing that the regular 17gon can be constructed by compass and ruler.
We implemented this algorithm of Gauss for computing radical expressions for roots of unity in MAPLE with an improved time complexity compared to Gauss' original proposition [Web96]. The command radsolve serves as an extension to the
solve
command, that returns radical solutions for all cyclotomic factors of a polynomial with integer coefficients and calls
solve
for the others.
The hard part of the proof and of the algorithm is to compute a radical expression for a primitive 
-th root of unity, where 
is a prime number. The improvement in our algorithm reduces the asymptotic time complexity in this case from
Ο(( 
to
Ο(
,
where is the largest prime factor of -1. An analysis of the algorithm and the statement and proof of the proposition on which the improvement is based can be found in [Web96]. We also refer to this paper for a description of the underlying ideas of Gauss' algorithm and the statement of the theorems that yield its correctness.
For an overview see The Algorithm.
Compared to the MAPLE-code given in [Web96], we managed to improve the practical speed of the algorithm to a great extend, reducing the needed amount of memory at the same time. For instance, the computations of radical expressions for the 47th, 59th, or 83rd root of unity failed with the previous implementation, but can now be accomplished within a few hours.
The Algorithm
First, we will show how the task of computing radical solutions for cyclotomic polynomials can be reduced to the one of calculating radical expressions for roots of unity. Second, we will sketch the algorithm for the latter task, which mainly is the one given in [Web96].
Radical solutions for cyclotomic polynomials
The n-th cyclotomic Polynomial over the rational numbers is of the form
![`𝓅`[n]](/view.aspx?SI=7067/radsolvelib_description_7.gif)
(x) = 
and its degree is
φ(n) = 
where n = ![product(`^`(`𝓅`[i], r[i]), i = 1 .. s)](/view.aspx?SI=7067/radsolvelib_description_10.gif)
and φ is the Eulerian φ -Function, as can be seen in [Nar90, p. 13 and 169]. From the difining formula of φ it can be seen that the inverse image 
of φ is finite. Thus, to determine whether or not an irreducible factor of a polynomial is cyclotomic, it suffices to compare the j-th cyclotomic polynomial with q for all j in 
. MAPLE provides the function numtheory[invphi] to calculate 
; the computational costs of this functions is neglegible for the range of arguments which are feasible for computing radical expressions.
Given the n-th cyclotomic polynomial ![`𝓅`[n]](/view.aspx?SI=7067/radsolvelib_description_14.gif)
(x) , the task of computing φ(n) radical solutions can be reduced to finding a radical expression for one root of unity, because this root will be primitive, e.g. a multiplicative generator for the others. Further, it suffices to find a primitive n-th root of unity to compute primitive roots of unity for all prime factors of n by applying some relatively simple recursion formulas, see e.g. [Gaa88].
Radical expressions for primitive r-th roots of unity
The hard part of the algorithm is to find a radical expression for a primitive 
-th root of unity, where 
is a prime number. This part involves the computation of Gaussian periods, cf. [Gau86, Sec.343]. If ![Zeta[`𝓅`]](/view.aspx?SI=7067/radsolvelib_description_17.gif)
denotes a primitive -th root of unity the period (f , k) for two integers f and k is defined by
(f , k) = 
where e = 
and ℊ is a primitive -th root. In Maple a primitive -th root can be computed by the function primroot in the package numtheory. The Gaussian periods are generators of intermediate fields of the number field Q( ![Zeta[`𝓅`]](/view.aspx?SI=7067/radsolvelib_description_20.gif)
), which is a number field of degree -1 over the rationals. Using the modern terminlogy of intermediate fields the algorithm of Gauss works roughly as follows. Starting with the rationals, compute a radical expression for a generator of an intermediate field of prime degree, i.e. an appropriate Gaussian period. After this task has been accomplished compute a radical expressions for a generator of another intermediate field (i.e. another Gaussian period) whose relative degree over the other intermediate field is of prime order. Continue this task until a radical expression for a generator of Q( ![Zeta[`𝓅`]](/view.aspx?SI=7067/radsolvelib_description_21.gif)
), i.e. ![Zeta[`𝓅`] = (1, `^`(`ℊ`, 0))](/view.aspx?SI=7067/radsolvelib_description_22.gif)
has been computed.
The difficult part of Gauss algorithm is the lifting of radical expressions from one intermediate field to another one. We will not give the details of this lifting here but we refer to [Web96].
The algorithm described above works for several sequences of intermediate fields to be constructed. All permutations on the list of prime factors of -1 give a possible sequence in the construction of the intermediate fields. Different permutations give different results in general, and also the computational costs are greatly affected by the choice of an appropriate sequence of intermediate fields. Our default strategy in the choice of the intermediate fields is as follows: first take the largest prime factor, then the second largest and so on, until we just get the final field Q( ![Zeta[`𝓅`]](/view.aspx?SI=7067/radsolvelib_description_23.gif)
) as a relative extension of degree 2 of the last intermediate field. Note that we have to count multiplicities as seperate factors. This strategy seems to be preferable over others, because bigger relative extensions --- which cause more computations --- occurr in smaller intermediate fields. However, as a tool for experimental mathematics we provide another strategy, too. The second strategy is simply to reverse the ordering; this can be done by setting SwapPrimeOrder := true. As had to be expected in all our experiments the obtained computation times were far worse than with our default setting and the computed radical expressions were also larger; some examples are given below.
How to use radsolve
The library is included in the file ?radsolvelib.mpl?.
All functionality is available via the function radsolve.
The command radsolve serves as an extension to the Maple solve command. It returns radical solutions for all cyclotomic factors of a polynomial with integer coefficients and calls solve for the others.
Description:
- Calling Sequence: radsolve(poly)
- Parameters: poly - an univariate polynomial with integer coefficients
- The command radsolve(poly) returns the solutions of the univariate polynomial poly as radical expressions for all the irreducible factors of poly that are either cyclotomic or of degree lower than four.
- In general, explicit solutions in terms of radicals for polynomial equations of degree greater than four do exist if and only if the Galois group of the polynomial is solvable. This is the case for cyclotomic polynomials, since their Galois group is cyclic. radsolve(poly) computes radical expressions for such cases and calls
solve
for the other factors of poly. In cases, where explicit solutions can't be computed, implicit solutions are given in terms of
RootOf
s.
- Note that to obtain explicit solutions for the general quartic polynomial equation you have to set the global variable _EnvExplicit to true, as explained in the topic help to
solve
.
- The output from radsolve is a sequence of solutions.
- The behaving of radsolve(poly) can be changed by using the following global variables (default values in brackets):
- _AllSolutions (false): If _AllSolutions is set to false, radsolve(poly) returns only one solution (no matter which exponent belongs to the factor) for every irreducible factor of poly . If it is set to true, all solutions are given.
- _UseRootsymbols (false): if _UseRootsymbols is set to true,
RootOf
s are used to represent primitive roots of unity of prime factors of , where n is the degree of the cyclotomic polynomial, although these roots could be expressed in radicals explicitly. This saves time and space. If _UseRootsymbols is set to false, all
RootOf
s are substituted by their radical expressions.
- _SwapPrimeOrder (false): If _SwapPrimeOrder is set to true, then an alternative order of building intermediate generators of field extensions of order , where is the degree of the cyclotomic polynomial and m divides , is chosen. This results in different radical expressions and in more computing time. Mainly for experimentation.
- To obtain information about the computing process, set infolevel[radsolve] to an appropriate level. Level 0 gives just the solutions, level 1 and level 2 are not used, level 3 gives a general outline and intermediate timing, level 4 is nice to use for larger examples and level 5 is designed for those who are really interested in the algorithm.
radsolve(poly) is based on an algorithm given by Carl Friedrich Gauss in his Disquisitiones Arithmeticae, 1797, with an improved time complexity, as described in Weber, A.: Computing radical expressions for roots of unity, SIGSAM Bulletin 30, 117 (Sept. 1996), p.11-20.
Examples:
The Maple solve command does not express the solutions of a cyclotomic factor of a polynomial of degree higher than 4 in radicals, but uses sin and cos functions instead.
In the following, we set AllSolutions to false. Thus we get only one solution for any irreducible factor of the polynomial.
| > |
fname1 := cat(libdir, "\\radsolvelib.mpl"):; |
| > |
; 1](/view.aspx?SI=7067/radsolvelib_description_35.gif) |
 |
(3.2.3) |
| > |
solve(16*x^3-16*x^6-14*x^2+2*x^9-12*x^8+13*x^7-3+4*x+16*x^5-16*x^4); |
| > |
 |
| > |
; 1](/view.aspx?SI=7067/radsolvelib_description_81.gif) |
 |
(3.2.6) |
| > |
infolevel[radsolve]:=3: |
| > |
radsolve(16*x^3-16*x^6-14*x^2+2*x^9-12*x^8+13*x^7-3+4*x+16*x^5-16*x^4); |
radsolve: The global variables are set to:
radsolve: _AllSolutions= false
radsolve: _SwapPrimeOrder= false
radsolve: _UseRootsymbols= true
radsolve: _EnvExplicit= false
radsolve: Degree less than 5 -> use solve
radsolve: invphi(degree)= [7 9 14 18]
radsolve: No, not equal.
radsolve: No, not equal.
radsolve: Yes, it's equal.
radsolve: It?s the 14 -th cyclotomic polynomial!
radsolve: Calculating the solution(s)...
radical_prime_primroot: Entering. Calculating a radical expression for a primitive 7 -th root of unity...
expressions_for_periods: Entering with prime factors of 6 : 3 2
expressions_for_periods: Entering with prime factors of 6 : 1 3
expressions_for_periods: Leaving. Generator for intermediate field of degree 3 computed.
expressions_for_periods: used time 0.15e-1
expressions_for_periods: Leaving. Generator for intermediate field of degree 6 computed.
expressions_for_periods: used time 0.31e-1
radical_prime_primroot: Entering. Calculating a radical expression for a primitive 3 -th root of unity...
expressions_for_periods: Entering with prime factors of 2 : 1 2
expressions_for_periods: Leaving. Generator for intermediate field of degree 2 computed.
expressions_for_periods: used time 0.
radical_prime_primroot: Calculated a radical expression for a primitive 3 -th root of unity. Leaving.
radical_prime_primroot: Calculated a radical expression for a primitive 7 -th root of unity. Leaving.
radical_primroot: Leaving.
radsolve: Used Time= 0.47e-1 |
 |
(3.2.7) |
See Also:
solve
,
RootOf
,
allvalues
,
dsolve
,
fsolve
,
isolve
,
msolve
,
rsolve
,
assign
,
invfunc
,
isolate
,
match
,
linalg[linsolve]
,
simplex
,
grobner
To verify the results, use
radnormal
. You might have to use numeric evaluation for large results:
| > |
z:=radsolve(numtheory[cyclotomic](41,x)): |
| > |
numeval:=codegen[makeproc](map(evalf,[codegen[optimize](z)])): |
| > |
fnormal(numeval()^41-1,30); |
radsolve: The global variables are set to:
radsolve: _AllSolutions= false |
radsolve: _SwapPrimeOrder= false
radsolve: _UseRootsymbols= true
radsolve: _EnvExplicit= false
radsolve: invphi(degree)= [41 55 75 82 88 100 110 132 150]
radsolve: Yes, it's equal.
radsolve: It?s the 41 -th cyclotomic polynomial!
radsolve: Calculating the solution(s)...
radical_prime_primroot: Entering. Calculating a radical expression for a primitive 41 -th root of unity...
expressions_for_periods: Entering with prime factors of 40 : 20 2
expressions_for_periods: Entering with prime factors of 40 : 10 2
expressions_for_periods: Entering with prime factors of 40 : 5 2
expressions_for_periods: Entering with prime factors of 40 : 1 5 |
expressions_for_periods: Leaving. Generator for intermediate field of degree 5 computed.
expressions_for_periods: used time .375 |
expressions_for_periods: Leaving. Generator for intermediate field of degree 10 computed.
expressions_for_periods: used time .484 |
expressions_for_periods: Leaving. Generator for intermediate field of degree 20 computed.
expressions_for_periods: used time .828 |
expressions_for_periods: Leaving. Generator for intermediate field of degree 40 computed.
expressions_for_periods: used time 1.203
radical_prime_primroot: Calculated a radical expression for a primitive 41 -th root of unity. Leaving.
radical_primroot: Leaving.
radsolve: Used Time= 1.203 |
 |
(3.3.1) |
Gauss computed the value of cos(2*pi/17) up to 30 digits in 1801, see [BH08]. We will compute the numeric value of the our radical expression for a 17th root of unity as a test.
| > |
z:=radsolve(numtheory[cyclotomic](17,x)): |
| > |
Digits:=32: infolevel[radsolve]:=0: |
| > |
numeval17:=codegen[makeproc](map(evalf,[codegen[optimize](z)])): |
| > |
 |
| radsolve: The global variables are set to: |
radsolve: _AllSolutions= false
radsolve: _SwapPrimeOrder= false
radsolve: _UseRootsymbols= true
radsolve: _EnvExplicit= false
radsolve: invphi(degree)= [17 32 34 40 48 60]
radsolve: Yes, it's equal.
radsolve: It?s the 17 -th cyclotomic polynomial!
radsolve: Calculating the solution(s)...
radsolve: Used Time= 0. |
 |
(3.3.2) |
| > |
 |
 |
(3.3.3) |
Note that radsolve uses remember tables. So, if you compute the same solution twice, the second will be taken from memory and changes to global variables will not affect the result. Try restart in this case.
Practical Limitations of the Algorithm
Compared to [Web96] the implementation of the algorithm has been optimized. For the table below we apply radsolve on all cyclotomic polynomials of (prime) degree up to 257. It can be seen that the computing time depends on the size of the largest prime factor of to a great extend, as had to be expected from the theoretical complexity analysis given in [Web96]. We also summarized sizes of normal and dag representation, f meaning rational, r radical operations and a the number of dag variables used. Moreover, we save the computed radical expressions in files with the name radprimroot###.m. We also test numerically that the computed expressions are roots of unity (of the right order).
| > |
restart: libdir := currentdir();
currentdir(libdir);
fname1 := cat(libdir, "\\radsolvelib.mpl"):
read (fname1):
_UseRootsymbols:=true:
_AllSolutions:=false:
_SwapPrimeOrder:=false:
Digits:=30;
additions:=`f`:
multiplications:=`f`:
divisions:=`f`:
functions:=`r`:
assignments:=`a`:
printf(`\n p p-1 time (in sec.) size of tree size of dag test`);
printf(`\n--------------------------------------------------------------------------------------------------`);
for i from 7 to 257 do
if isprime(i) then
t:=round(time(radsolve(numtheory[cyclotomic](i,x)))):
z:= radsolve(numtheory[cyclotomic](i,x));
rname:=radprimroot||i;
rpr||i:=z;
save(rpr||i,cat(rname,".m"));
numeval:=codegen[makeproc](map(evalf,[codegen[optimize](z)])):
printf(`\n%5d %15A %10d %20A %25A %f`,i, ifactor(i-1),t,codegen[cost](z),codegen[cost](codegen[optimize](z)), evalf(Re(evalf(numeval()^i,25)),18));
fi:
od:
printf(`\n--------------------------------------------------------------------------------------------------`);
|
 |
 |
 |
p p-1 time (in sec.) size of tree size of dag test
--------------------------------------------------------------------------------------------------
7 (2)*(3) 0 95*f+24*r 5*r+28*f+7*a 1.000000
11 (2)*(5) 0 869*f+78*r 5*r+75*f+15*a 1.000000
13 (2)^2*(3) 0 291*f+74*r 9*r+52*f+13*a 1.000000
17 (2)^4 0 55*f+38*r 12*r+31*f+10*a 1.000000
19 (2)*(3)^2 0 1127*f+264*r 9*r+89*f+24*a 1.000000
23 (2)*(11) 4 25541*f+442*r 5*r+361*f+48*a 1.000000
29 (2)^2*(7) 1 11043*f+498*r 9*r+222*f+30*a 1.000000
31 (2)*(3)*(5) 1 12697*f+1138*r 10*r+201*f+36*a 1.000000
37 (2)^2*(3)^2 1 3567*f+822*r 15*r+163*f+29*a 1.000000
41 (2)^3*(5) 1 10677*f+958*r 19*r+234*f+75*a 1.000000
43 (2)*(3)*(7) 2 59778*f+2688*r 10*r+354*f+59*a 1.000000
47 (2)*(23) 71 515406*f+2026*r 5*r+1409*f+113*a 1.000000
53 (2)^2*(13) 14 151629*f+1878*r 9*r+656*f+85*a 1.000000
59 (2)*(29) 195 1322399*f+3250*r 5*r+2192*f+136*a 1.000000
61 (2)^2*(3)*(5) 3 41184*f+3730*r 20*r+337*f+66*a 1.000000
67 (2)*(3)*(11) 12 400677*f+7052*r 10*r+695*f+129*a 1.000000
71 (2)*(5)*(7) 8 185432*f+8374*r 10*r+599*f+101*a 1.000000
73 (2)^3*(3)^2 3 17788*f+4052*r 25*r+289*f+77*a 1.000000
79 (2)*(3)*(13) 23 798813*f+10002*r 10*r+917*f+136*a 1.000000
83 (2)*(41) 827 5393221*f+6562*r 5*r+4399*f+208*a 1.000000
89 (2)^3*(11) 16 376715*f+6608*r 19*r+893*f+203*a 1.000000
97 (2)^5*(3) 6 10015*f+2566*r 37*r+222*f+55*a 1.000000
101 (2)^2*(5)^2 11 106153*f+9194*r 19*r+685*f+85*a 1.000000
103 (2)*(3)*(17) 58 2408025*f+17402*r 10*r+1378*f+180*a 1.000000
107 (2)*(53) 2594 15333111*f+11026*r 5*r+7357*f+375*a 1.000000
109 (2)^2*(3)^3 8 94540*f+21262*r 21*r+403*f+84*a 1.000000
113 (2)^4*(7) 11 270760*f+12042*r 33*r+742*f+195*a 1.000000
127 (2)*(3)^2*(7) 16 1268882*f+57204*r 24*r+829*f+164*a 1.000000
131 (2)*(5)*(13) 48 4205545*f+52882*r 10*r+1408*f+252*a 1.000000
137 (2)^3*(17) 79 2268411*f+16338*r 19*r+1740*f+300*a 1.000000
139 (2)*(3)*(23) 198 8276511*f+32382*r 10*r+2192*f+284*a 1.000000
149 (2)^2*(37) 961 10728147*f+15954*r 9*r+4109*f+278*a 1.000000
151 (2)*(3)*(5)^2 25 584637*f+50574*r 18*r+848*f+132*a 1.000000
157 (2)^2*(3)*(13) 51 4952625*f+61876*r 20*r+1654*f+301*a 1.000000
163 (2)*(3)^4 26 527679*f+118672*r 31*r+746*f+139*a 1.000000
167 (2)*(83) 20402 92925886*f+27226*r 5*r+17801*f+586*a 1.000000
173 (2)^2*(43) 1999 19641408*f+21678*r 9*r+5387*f+347*a 1.000000
179 (2)*(89) 24954 123295437*f+31330*r 5*r+20441*f+979*a 1.000000
181 (2)^2*(3)^2*(5) 30 1616894*f+146440*r 30*r+1017*f+166*a 1.000000
191 (2)*(5)*(19) 179 27351279*f+157930*r 10*r+2419*f+353*a 1.000000
193 (2)^6*(3) 29 103141*f+25926*r 61*r+485*f+112*a 1.000000
197 (2)^2*(7)^2 57 955995*f+41994*r 19*r+1676*f+203*a 1.000000
199 (2)*(3)^2*(11) 59 12001853*f+211954*r 26*r+1567*f+351*a 1.000000
211 (2)*(3)*(5)*(7) 57 3089840*f+139716*r 23*r+1396*f+194*a 1.000000
223 (2)*(3)*(37) 1345 57178095*f+85290*r 10*r+4971*f+451*a 1.000000
227 (2)*(113) 73657 321494265*f+50626*r 5*r+32736*f+1250*a 1.000000
229 (2)^2*(3)*(19) 200 30660759*f+176446*r 20*r+3164*f+516*a 1.000000
233 (2)^3*(29) 656 19989027*f+48710*r 19*r+3893*f+525*a 1.000000
239 (2)*(7)*(17) 278 25676822*f+186676*r 10*r+2838*f+335*a 1.000000
241 (2)^4*(3)*(5) 59 1914333*f+173866*r 48*r+883*f+223*a 1.000000
251 (2)*(5)^3 129 5483003*f+473468*r 19*r+1697*f+272*a 1.000000
257 (2)^8 60 12333*f+8390*r 56*r+289*f+84*a 1.000000
-------------------------------------------------------------------------------------------------- |
|
It is interesting to see the results you obtain by setting SwapPrimeOrder := true. You get different radical expressions and it takes more time and memory. The reason for this is explained above. We will only give the dag representation. Notice that you must not use a simply evalf for the numerical test, as this would generate a tree representation of the expression instead of a dag.
| > |
restart: libdir := currentdir();
currentdir(libdir);
fname1 := cat(libdir, "\\radsolvelib.mpl"):
read (fname1):
_UseRootsymbols:=true:
_AllSolutions:=false:
_SwapPrimeOrder:=true:
Digits:=30;
additions:=`f`:
multiplications:=`f`:
divisions:=`f`:
functions:=`r`:
assignments:=`a`:
printf(`\n p p-1 time (in sec.) size of dag test`);
printf(`\n--------------------------------------------------------------------------------------------------`);
plist:=[23,43,71]:
for i in plist do
t:=round(time(radsolve(numtheory[cyclotomic](i,x)))):
z:= radsolve(numtheory[cyclotomic](i,x));
rname:=radprimrootalt||i;rpr||i:=z;
save(rpr||i,cat(rname,".m"));
numeval:=codegen[makeproc](map(evalf,[codegen[optimize](z)])):
printf(`\n%5d %15A %10d %25A %f`,i, ifactor(i-1),t,codegen[cost](codegen[optimize](z)), evalf(Re(evalf(numeval()^i,25)),18));
od: |
We will also compute the fully recursively evaluated radical expressions for the prime roots of unity. Here we will only give the count on the dag representation. Notice that we must not use evalf directly, as this will cause a possible explosion on memory due to an intermediate tree representation of the expression. We save the results in files of the name radicalprimroot###.m.
 |
 |
 |
p p-1 time (in sec.) size of dag test
--------------------------------------------------------------------------------------------------
23 (2)*(11) 6 5*r+723*f+48*a 1.000000
43 (2)*(3)*(7) 7 10*r+547*f+55*a 1.000000
71 (2)*(5)*(7) 29 10*r+938*f+106*a 1.000000 |
|
| > |
restart: libdir := currentdir();
currentdir(libdir); fname1 := cat(libdir, "\\radsolvelib.mpl"):
read (fname1):
_UseRootsymbols:=false:
_AllSolutions:=false:
_SwapPrimeOrder:=false:
additions:=`f`:
multiplications:=`f`:
divisions:=`f`:
functions:=`r`:
assignments:=`a`:Digits:=30;
printf(`\n p p-1 time (in sec.) size of dag test`);
printf(`\n--------------------------------------------------------------------------------------------------`);
for i from 7 to 257 do
if isprime(i) then
t:=round(time(radsolve(numtheory[cyclotomic](i,x)))):
z:= radsolve(numtheory[cyclotomic](i,x));
rname:=radicalprimroot||i;rprr||i:=z;
save(rprr||i,cat(rname,".m"));
numeval:=codegen[makeproc](map(evalf,[codegen[optimize](z)])):
printf(`\n%5d %15A %10d %30A %f`,i, ifactor(i-1),t,codegen[cost](codegen[optimize](z)), evalf(Re(evalf(numeval()^i,25)),18));
fi:
od:
printf(`\n--------------------------------------------------------------------------------------------------`); |
 |
 |
 |
p p-1 time (in sec.) size of dag test
--------------------------------------------------------------------------------------------------
7 (2)*(3) 0 6*r+27*f+7*a 1.000000
11 (2)*(5) 0 8*r+79*f+18*a 1.000000
13 (2)^2*(3) 0 10*r+53*f+14*a 1.000000
17 (2)^4 0 12*r+31*f+10*a 1.000000
19 (2)*(3)^2 0 10*r+90*f+24*a 1.000000
23 (2)*(11) 3 12*r+460*f+72*a 1.000000
29 (2)^2*(7) 1 14*r+248*f+40*a 1.000000
31 (2)*(3)*(5) 1 14*r+211*f+40*a 1.000000
37 (2)^2*(3)^2 1 16*r+164*f+31*a 1.000000
41 (2)^3*(5) 1 22*r+239*f+78*a 1.000000
43 (2)*(3)*(7) 2 14*r+381*f+65*a 1.000000
47 (2)*(23) 89 16*r+2060*f+199*a 1.000000
53 (2)^2*(13) 13 18*r+702*f+98*a 1.000000
59 (2)*(29) 182 18*r+2423*f+194*a 1.000000
61 (2)^2*(3)*(5) 3 24*r+342*f+71*a 1.000000
67 (2)*(3)*(11) 11 18*r+804*f+152*a 1.000000
71 (2)*(5)*(7) 7 18*r+665*f+114*a 1.000000
73 (2)^3*(3)^2 3 26*r+294*f+78*a 1.000000
79 (2)*(3)*(13) 21 18*r+959*f+151*a 1.000000
83 (2)*(41) 799 26*r+4284*f+367*a 1.000000
89 (2)^3*(11) 16 26*r+1016*f+230*a 1.000000
97 (2)^5*(3) 6 38*r+223*f+56*a 1.000000
101 (2)^2*(5)^2 11 22*r+739*f+102*a 1.000000
103 (2)*(3)*(17) 56 22*r+1374*f+191*a 1.000000
107 (2)*(53) 2735 22*r+8034*f+486*a 1.000000
109 (2)^2*(3)^3 8 22*r+401*f+89*a 1.000000
113 (2)^4*(7) 11 38*r+770*f+203*a 1.000000
127 (2)*(3)^2*(7) 16 28*r+852*f+172*a 1.000000
131 (2)*(5)*(13) 47 22*r+1497*f+278*a 1.000000
137 (2)^3*(17) 77 30*r+1754*f+330*a 1.000000
139 (2)*(3)*(23) 529 22*r+2850*f+380*a 1.000000
149 (2)^2*(37) 942 24*r+3988*f+313*a 1.000000
151 (2)*(3)*(5)^2 24 22*r+903*f+149*a 1.000000
157 (2)^2*(3)*(13) 50 28*r+1701*f+312*a 1.000000
163 (2)*(3)^4 25 32*r+749*f+143*a 1.000000
167 (2)*(83) 77555 30*r+23002*f+1273*a 1.000000
173 (2)^2*(43) 2148 22*r+5538*f+460*a 1.000000
179 (2)*(89) 29892 30*r+19652*f+1208*a 1.000000
181 (2)^2*(3)^2*(5) 30 34*r+1044*f+175*a 1.000000
191 (2)*(5)*(19) 180 22*r+2543*f+404*a 1.000000
193 (2)^6*(3) 28 62*r+486*f+113*a 1.000000
197 (2)^2*(7)^2 57 24*r+1899*f+222*a 1.000000
199 (2)*(3)^2*(11) 62 34*r+1677*f+376*a 1.000000
211 (2)*(3)*(5)*(7) 56 30*r+1491*f+210*a 1.000000
223 (2)*(3)*(37) 1332 24*r+4815*f+477*a 1.000000
227 (2)*(113) 75056 42*r+29182*f+1466*a 1.000000
229 (2)^2*(3)*(19) 206 28*r+3194*f+536*a 1.000000
233 (2)^3*(29) 716 32*r+4316*f+575*a 1.000000
239 (2)*(7)*(17) 277 26*r+3131*f+370*a 1.000000
241 (2)^4*(3)*(5) 59 52*r+881*f+226*a 1.000000
251 (2)*(5)^3 129 22*r+1813*f+296*a 1.000000
257 (2)^8 60 56*r+289*f+84*a 1.000000
-------------------------------------------------------------------------------------------------- |
|
Comparison to Related Work
In the book of Gaal [Gaa88] a radical expression for a 7th root of unity is derived by some special reasoning that does not generalize to higher orders. The derived expression is the following:
 |
 |
(5.1) |
| > |
A:=-1/6+(7/2+21/2*sqrt(-3))^(1/3)/6+(7/2-21/2*sqrt(-3))^(1/3)/6+(sqrt((-1/3+(7/2+21/2*sqrt(-3))^(1/3)/3+(7/2-21/2*sqrt(-3))^(1/3)/3)^2-4))/2; |
 |
(5.2) |
Our function radsolve computes:
| > |
B:=radsolve(x^6+x^5+x^4+x^3+x^2+x+1); |
radsolve computed the 7th root of unity that is equal to 
:
 |
(5.4) |
Except for an implementation by ourselves of a much more inefficient method for the same task described in [Edw84], we do not know of implementations of other general methods by which radical expressions for higher roots of unity can be computed. Without being aware of an implemetation, we know of an algorithm developped by B. Trager, which computes radical expressions for a -th root of unity [Zip94]. This algorithm is entirely different from the one of Gauss. The major computational task consists of inverting a matrix of size Ο() over Q( ![Zeta[`𝓅`]](/view.aspx?SI=7067/radsolvelib_description_119.gif)
) , where q is a divisor of -1. Thus if -1 is smooth, i.e. if -1 contains only small prime factors, the asymptotic time complexity of our improvement of the algorithm of Gauss is much better. But in special cases, such that 
is prime, the algorithm of Trager might be an interesting alternative. It would be interesting to have a careful implementation of the hard cases, such as p = 47 , 59, 83, 107, 167, 179, 227 etc.
Bibliography
[BH08] Bauer, F. L. and Haenel, C.: Carl Friedrich Gau?, das 17-Eck und MATHEMATICA. Informatik Spektrum 31,5 (2008), p. 492-498
[Edw84] Edwards, H. M.: Galois Theory, vol. 101 of Graduate Texts in Mathematics, Springer, New York, 1984
[Gaa88] Gaal, L.: Classical Galois Theory, 4th ed., Chelsea Publishing Company, New York, 1988
[Gau86] Gauss, C. F.: Disquisitiones Arithmeticae - English Edition, Springer, Berlin, 1986
[HM01] Hanrot, G. and Morain, F.: Solvability of radicals from an algorithmic point of view. Proc. ISSAC 2001, p. 175-182. ACM.
[Nar90] Narkiewicz, W.: Elementary and Analytic Theory of Numbers, Springer, 1990.
[Web96] Weber, A.: Computing radical expressions for roots of unity, SIGSAM Bulletin 30, 117 (Sept. 1996), p. 11-20
[Zip94] Zippel, R.: Computer Algebra. Unpublished Lecture Notes, 1994.
