I only now received this enormously helpful reply that Julian sent 11 days ago! -------------------------- This works by guessing that the quintic has solutions of the form c+w*r1^(1/5)+w^2*r2^(1/5)+w^3*r3^(1/5)+w^4*r4^(1/5), where r1, r2, r3, r4 are the roots of a quartic and w ranges over fifth roots of unity. The inverse Fourier transform is how you go from the set of roots to (the fifth roots of) r1, r2, r3, r4. The rest of the code is just trying to guess the right order to put the roots in and which fifth roots of unity to use. This works if and only if the solutions have the form I described (and the numerics are good enough), and algebraists know (but I don't) whether that covers all solvable cases—it's a fairly simple Galois theory question. Julian ------------------------ Since squint[q_] (defined below) is so enormously simpler than published (Galois-theoretic) algorithms, one has to suspect that it is either unknown, or known to fail. If it is merely *conjectured* to fail, that would seem worthy of investigation. Note that it doesn't even bother with the usual solvability checks. So far, it has failed to fail on hundreds of cases. Here's what happens with 2x⁵ -5x-4. (Can anybody guess why I remember this as JRST?) First, numerically approximate the roots: In[649]:= NSolve[4 == 2 x^5 - 5 x] Out[649]= {{x -> -0.871327718648117 - 0.225067056441812 I}, {x -> -0.871327718648117 + 0.225067056441812 I}, {x -> 0.167726766484069 - 1.31407674104546 I}, {x -> 0.167726766484069 + 1.31407674104546 I}, {x -> 1.4072019043281}} We need to search over the 5! = 120 permutations for one of the 20 which actually works. I peeked ahead and found that the first two are the 2nd and 9th in Mathematica's preferred order: In[674]:= Permutations[%651][[{2, 9}]] Out[674]= {{-0.871327718648117 - 0.225067056441812 I, -0.871327718648117 + 0.225067056441812 I, 0.167726766484069 - 1.31407674104546 I, 1.4072019043281, 0.167726766484069 + 1.31407674104546 I}, {-0.871327718648117 - 0.225067056441812 I, 0.167726766484069 - 1.31407674104546 I, 0.167726766484069 + 1.31407674104546 I, -0.871327718648117 + 0.225067056441812 I, 1.4072019043281}} In[678]:= inversefourier /@ %674/√5 Out[678]= {{0. + 0. I, -0.8342002230717 + 0.606081938933149 I, 0.158846909677956 - 0.48888051884701 I, -0.0850160564298369 - 0.261652517269671 I, -0.110958348824536 - 0.0806159592582794 I}, {0. - 1.66533453693773*10^-17 I, -0.0850160564298369 - 0.261652517269671 I, -0.8342002230717 + 0.606081938933149 I, -0.110958348824536 - 0.0806159592582794 I, 0.158846909677956 - 0.48888051884701 I}} In[679]:= %^5; In[680]:= Chop@% Out[680]= {{0, 1.16563689479139, 0.0358907084015456, -0.00157613335078357, 0.0000485301578456363}, {0, -0.00157613335078357, 1.16563689479139, 0.0000485301578456363, 0.0358907084015455}} Permute[Fourier] = Fourier[Permute] ?? Maybe that's a criterion. Form the symmetric polynomials: In[681]:= Times @@ (x - Rest@#) & /@ % Out[681]= {(-1.16563689479139 + x) (-0.0358907084015456 + x) (-0.0000485301578456363 + x) (0.00157613335078357 + x), (-1.16563689479139 + x) (-0.0358907084015455 + x) (-0.0000485301578456363 + x) (0.00157613335078357 + x)} In[682]:= Rationalize@Expand@% Out[682]= {-(1/312500000) + x/15625 + x^2/25 - (6 x^3)/5 + x^4, -(1/312500000) + x/ 15625 + x^2/25 - (6 x^3)/5 + x^4} A quartic with apparently rational coefficients. Twice! Start getting excited. If the original polynomial 2x⁵ - 5x -4 was ax⁵+bx⁴+..., compute b/a: In[683]:= boa = Rationalize@Mean[x /. %649] Out[683]= 0 Solve the quartic and take 5th roots: In[684]:= Together[#1[[1, 2]]]^(1/5) & /@ Solve[%682[[1]] == 0] Out[684]= {Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 1, 0]^(1/5), Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 2, 0]^(1/5), Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 3, 0]^(1/5), Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 4, 0]^(1/5)} Then the grand solution is the union of all these tuples that work numerically: In[685]:= boa + Union[(% . #1 & ) /@ Select[Tuples[I^(4 Range@5/5), 4], MemberQ[Chop[#1 . % + boa - x /. %649], 0] & ]] Out[685]= {(-1)^(2/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 1, 0]^(1/5) + (-1)^(4/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 2, 0]^(1/5) + (-1)^(2/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 3, 0]^(1/5) - \ (-Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 4, 0])^( 1/5), (-1)^(4/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 1, 0]^(1/5) + Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 2, 0]^( 1/5) + Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 3, 0]^(1/5) + Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 4, 0]^( 1/5), Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 1, 0]^(1/5) - (-1)^(3/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 2, 0]^(1/5) + (-1)^(4/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 3, 0]^ (1/5) + (-1)^(2/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 4, 0]^(1/5), -(-1)^(3/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 1, 0]^(1/5) + (-1)^(2/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 2, 0]^(1/5) - \ (-Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 3, 0])^( 1/5) - (-1)^(3/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 4, 0]^(1/5), -(-1)^(1/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 1, 0]^(1/5) - \ (-Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 2, 0])^( 1/5) - (-1)^(3/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 3, 0]^(1/5) + (-1)^(4/5) Root[-1 + 20000 # + 12500000 #^2 - 375000000 #^3 + 312500000 #^4& , 4, 0]^(1/5)} Test the first one: In[688]:= RootReduce@%685[[1]] // tim During evaluation of In[688]:= 23.119116,3 Out[688]= Root[-4 - 5 # + 2 #^5& , 3, 0] Bingo! Grab and test the rest: In[691]:= tim@MinimalPolynomial@# & /@ %685 During evaluation of In[691]:= 12.150844,1 (seconds) During evaluation of In[691]:= 15.811916,1 During evaluation of In[691]:= 15.652511,1 During evaluation of In[691]:= 17.422262,1 During evaluation of In[691]:= 13.513089,1 Out[691]= {-4 - 5 #1 + 2 #1^5 &, -4 - 5 #1 + 2 #1^5 &, -4 - 5 #1 + 2 #1^5 &, -4 - 5 #1 + 2 #1^5 &, -4 - 5 #1 + 2 #1^5 &} It seems a little bassackwards to check by reconstructing the polynomial from the radicals, instead of plugging the radicals into the polynomial, but I forget how to goad Mathematica into this in finite time. You might object that the Root expressions in Out[685] are not in radicals. But squiint's outerrmost function is ToRadicals, which will indeed "radicalize" the quartics. In[737]:= #@%648 & /@ {LeafCount, ByteCount} Out[737]= {857, 25768} I.e., the full solution is almost 26KB. Little post-simplification is possible. squint generally returns short answers when possible. Credit to Mathematica's quartic solver. As always, Julian is right about the numerics. This simplistic implementation will eventually fail on problems with huge coefficients. —rwg On Mon, Nov 2, 2020 at 6:51 AM Bill Gosper <billgosper@gmail.com> wrote:
Unfortunately, I no longer have enough brain cells to remember or figure out how it works. It somehow gets exact expressions by taking the inverse Fourier transform of numerically approximated roots, without all the Galois group and resolvent rigamarole. I.e., it claims to find the radicals when they exist, subject to adjustable limits on coefficient size. Two examples (tim is just a timing function):
In[451]:= squint[x^5 - 5 x x - 3] // tim
During evaluation of In[451]:= 0.319225 seconds, 5 roots
Out[451]= {-((-1)^(3/5)/(2^(2/5) (5/(5 - Sqrt[5] - Sqrt[6 (5 - 11/Sqrt[5])]))^( 1/5))) + ((-1)^(2/5) (1/5 (5 - Sqrt[5] + Sqrt[6 (5 - 11/Sqrt[5])]))^(1/5))/2^( 2/5) - (-1 - 1/Sqrt[5] - 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5)/2^( 2/5) - ((-1)^(3/5) (1 + 1/Sqrt[5] - 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5))/2^(2/5),
1/(2^(2/5) (5/(5 - Sqrt[5] - Sqrt[6 (5 - 11/Sqrt[5])]))^( 1/5)) + (1/5 (5 - Sqrt[5] + Sqrt[6 (5 - 11/Sqrt[5])]))^(1/5)/2^( 2/5) + ((-1)^(4/5) (1 + 1/Sqrt[5] - 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5))/2^( 2/5) + (1 + 1/Sqrt[5] + 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5)/2^(2/5),
((-1)^(4/5) (1/5 (5 - Sqrt[5] + Sqrt[6 (5 - 11/Sqrt[5])]))^( 1/5))/2^(2/5) - (1/5 (-5 + Sqrt[5] + Sqrt[6 (5 - 11/Sqrt[5])]))^( 1/5)/2^(2/5) + ((-1)^(2/5) (1 + 1/Sqrt[5] - 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5))/2^( 2/5) + ((-1)^(2/5) (1 + 1/Sqrt[5] + 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5))/2^(2/5),
(-1)^(4/5)/(2^(2/5) (5/(5 - Sqrt[5] - Sqrt[6 (5 - 11/Sqrt[5])]))^( 1/5)) - (1/5 (-5 + Sqrt[5] - Sqrt[6 (5 - 11/Sqrt[5])]))^(1/5)/2^( 2/5) - ((-1)^(1/5) (1 + 1/Sqrt[5] - 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5))/2^( 2/5) - ((-1)^(3/5) (1 + 1/Sqrt[5] + 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5))/2^(2/5),
(-1)^(2/5)/( 2^(2/5) (5/(5 - Sqrt[5] - Sqrt[6 (5 - 11/Sqrt[5])]))^( 1/5)) - ((-1)^(3/5) (1/5 (5 - Sqrt[5] + Sqrt[6 (5 - 11/Sqrt[5])]))^(1/5))/2^( 2/5) + (1 + 1/Sqrt[5] - 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5)/2^( 2/5) + ((-1)^(4/5) (1 + 1/Sqrt[5] + 2 Sqrt[3/10 + 33/(50 Sqrt[5])])^(1/5))/2^(2/5)}
In[445]:= squint[-1 + 13755 x + 3740 x^2 + 1260 x^3 - 5 x^4 + x^5]
Out[445]= {1 + 5 2^(1/5) - 10 2^(2/5) + 10 2^(3/5) - 5 2^(4/5),
1 + 5 (-1)^(4/5) 2^(1/5) + 10 (-1)^(3/5) 2^(2/5) + 10 (-1)^(2/5) 2^(3/5) + 5 (-1)^(1/5) 2^(4/5),
1 - 5 (-1)^(3/5) 2^(1/5) + 10 (-1)^(1/5) 2^(2/5) + 10 (-1)^(4/5) 2^(3/5) - 5 (-1)^(2/5) 2^(4/5),
1 + 5 (-1)^(2/5) 2^(1/5) - 10 (-1)^(4/5) 2^(2/5) - 10 (-1)^(1/5) 2^(3/5) + 5 (-1)^(3/5) 2^(4/5),
1 - 5 (-2)^(1/5) - 10 (-2)^(3/5) - 10 (-1)^(2/5) 2^(2/5) - 5 (-1)^(4/5) 2^(4/5)}
(It can take Mathematica several minutes to check one of these symbolically.)
Here is the definition:
squint[q_] := ToRadicals[ Block[{rts = (#1[[1, 2]] & ) /@ NSolve[q == 0, WorkingPrecision -> 99], boa, foos}, foos = (InverseFourier[rts[[#1]]/Sqrt[5]]^5 & ) /@ Permutations[Range[5]]; foos = (Together[#1[[1, 2]]]^(1/5) & ) /@ Solve[0 == Rationalize[Expand[Times @@ (#1 - Rest[MinimalBy[foos, Denominator@Rationalize[Plus @@ #, 9.^-69] &,1][[1]]])]]]; boa = Rationalize@Mean@rts; boa + Union[(foos . #1 & ) /@ Select[Tuples[I^(4 Range@5/5), 4], MemberQ[Chop[#1 . foos + boa - rts], 0] & ]]]]
One thing I remember is boa is "b over a", as in a x^5 + b x^4 + . . . note boa = Rationalize@Mean@rts; Something I don't remember is feeling any surprise that this works. It would be funny if it comes as a surprise to algebraists. Actually, it wouldn't be funny because then I'd have to publish something. That I don't remember.
Replacing 5 by 7, 4 by 6, and Solve[0==. . .] with a sextic solver, Julian and I found the analogous septic solver. He expressed mistrust, but it has yet to fail. —rwg