冪乗に分解できるタイプ

今回のオリジナルプログラムでは，例えば，SolveSolvable(x^2-2)$ の根の表示には [[alpha[1],alpha[1]^2-8]] が現れるので，SolveSolvable(x^2-8)$ と問うと，[[alpha[1],alpha[1]^2-32]] が現れるので，．．．となってしまいます．

プログラムの趣意に反する可能性もありますが，例えば，x^2-2,x^6+x^3+1,x^20+x^15+3*x^5+4 といった冪乗に分解できるタイプについては，それぞれ x-2,x^2+x+1,x^4+x^3+3*x+4 の根の表示を得た上で，その各根の2,3,5乗根のすべてを出力するコードを書いてみました．

まず，1の原始N乗根の最小多項式を出力する関数を用意します．

ROU(x,N):=block([D:reverse(listify(divisors(N)))],
 num(factor((x^pop(D)-1)/apply("*",map(lambda([s],x^s-1),D)))))$
/* テスト */
makelist(ROU(x,i),i,1,10);

次に SolveSolvable の定義の局所変数の宣言をコメントアウトした SolveSolvable2.mac をリロードして，以下を実行すると，上記の分解を挟んだ結果が得られます．

/* サンプル選択 */
p:[x^2-2,x^6+x^3+1,x^20+x^15+3*x^5+4][2];
/* 分解（polydecompの出力は気まぐれなので．．．） */
for i:(N:hipow(p,x)) step -1 unless polynomialp(q:rat(subst(x=x^(1/i),p)),[x]) do N:i-1$[q,x^N];
/* ベースの方程式の求根 */
if hipow(q,x)=1 then [C,SolN,x[1]]:[[],1,rhs(solve(q,x)[1])] else (SolveSolvable(q),C:SI[StageN]@ExtensionList)$
/* 各根x[i]のN乗根y[i]に1の原始N乗根の冪乗を掛ける（得られるもの全体はy[i]の選び方に拠らない） */
for i:1 thru SolN do C:push([y[i],y[i]^N-x[i]],C)$push([w,ROU(w,N)],C);
RROOTS:flatten(makelist(makelist(y[i]*w^j,j,0,N-1),i,1,SolN));
/* 検算 */
ef_polynomial_reduction(P:apply("*",map(lambda([s],x-s),RROOTS)),C);