７．円分多項式、円分数

円分多項式

円分多項式Φ_n(x) とは、以下の定義により帰納的に定義される多項式である。

• ｎ＝１のとき、Φ₁(x)＝ｘ−１
• ｎ≧２のとき、

ｘn−１＝

Π d｜n

Φd(x)

（ｎの全ての約数ｄに関する積）。

円分多項式は、Z[X]上既約である。また、次数は、φ(n)（φはオイラー関数）となる。

Φ_n(x) については、以下のような直接的に表す式がある。

Φn(x)＝

Π d｜n

（ｘd−１）μ(n/d)

μ(n/d) はメビウス関数を表す。

この式から、円分多項式Φ_n(x) を求めてみる。

まず、ｎを割り切るような全てのｄによる積については、今、ｎの値が小さいので、
１≦ｄ≦ｎの範囲について、ｎを割り切るかどうかを直接調べていっても、それほど時間はかからない。

積については、メビウス関数の値に応じて、ｘ^d−１を掛ける場合と、ｘ^d−１で割る場合がある。
これについては、掛ける方を先に行う（途中の状態では割り切れるとは限らないので）。

プログラムは以下のとおり。

 10   'poly.ub
 20   L%=1000:dim A%(L%),B%(L%),W%(L%)
 30   input N%
 40   '
 50   ' main loop
 60   '
 70   Ka%=0:A%(0)=1
 80   for D%=1 to N%
 90     if or{N%@D%>0,moeb(N%\D%)<=0} then 110
100     gosub 210:gosub 250
110   next D%
120   for D%=1 to N%
130     if or{N%@D%>0,moeb(N%\D%)>=0} then 150
140     gosub 210:gosub 360
150   next D%
160   gosub 480:' print
170   print:goto 30
180   '
190   ' set (x^d - 1)
200   '
210   Kb%=D%:B%(Kb%)=1:B%(0)=-1:for I%=1 to Kb%-1:B%(I%)=0:next I%:return
220   '
230   ' w = a * b
240   '
250   Kw%=Ka%+Kb%:for I%=0 to Kw%:W%(I%)=0:next I%
260   for I%=0 to Ka%
270     for J%=0 to Kb%:K%=I%+J%
280       W%(K%)=W%(K%)+A%(I%)*B%(J%)
290     next J%
300   next I%
310   Ka%=Kw%:for I%=0 to Ka%:A%(I%)=W%(I%):next I%:return
320   return
330   '
340   ' w = a / b
350   '
360   Kw%=Ka%-Kb%
370   for I%=Kw% to 0 step -1
380     W%(I%)=A%(I%+Kb%)\B%(Kb%)
390     for J%=I% to I%+Kb%
400       A%(J%)=A%(J%)-W%(I%)*B%(J%-I%)
410     next J%
420   next I%
430   Ka%=Kw%:for I%=0 to Ka%:A%(I%)=W%(I%):next I%:return
440   return

プログラムについて簡単に説明すると、

• 80行目〜150行目がメインループ。
  80行目〜110行目が、メビウス関数の値が＋の場合、
  120行目〜150行目が、−の場合、となっている。
• 250行目〜320行目が多項式の掛算、360行目〜440行目が多項式の割算のルーチンである。

これで、計算が終わったところで a%(i%) に円分多項式の係数がセットされるので、
そのまま印字してもよいが、少し出力に凝ってみる。

450   '
460   ' pretty print
470   '
480   for I%=Ka% to 0 step -1
490     if A%(I%)<>0 then cancel for:goto 510
500   next I%
510   if I%=0 then print A%(I%);:return
520   if A%(I%)<>1 then print A%(I%);
530   print "x";
540   if I%>1 then print "^";cutspc(str(I%));
550   for J%=I%-1 to 0 step -1
560   ' coefficient
570     if A%(J%)=0 then 660
580     if A%(J%)>1 then print "+";cutspc(str(A%(J%)));
590     if A%(J%)=1 then print "+";
600     if A%(J%)=-1 then print "-";
610     if A%(J%)<-1 then print "-";cutspc(str(abs(A%(J%))));
620     if J%=0 then print cutspc(str(abs(A%(0))));
630   ' term
640     if J%>0 then print "x";
650     if J%>1 then print "^";cutspc(str(J%));
660   next J%:return

ｎ＝１〜20についての円分多項式は以下のとおり。

 1 : x-1
 2 : x+1
 3 : x^2+x+1
 4 : x^2+1
 5 : x^4+x^3+x^2+x+1
 6 : x^2-x+1
 7 : x^6+x^5+x^4+x^3+x^2+x+1
 8 : x^4+1
 9 : x^6+x^3+1
10 : x^4-x^3+x^2-x+1
11 : x^10+x^9+x^8+x^7+x^6+x^5+x^4+x^3+x^2+x+1
12 : x^4-x^2+1
13 : x^12+x^11+x^10+x^9+x^8+x^7+x^6+x^5+x^4+x^3+x^2+x+1
14 : x^6-x^5+x^4-x^3+x^2-x+1
15 : x^8-x^7+x^5-x^4+x^3-x+1
16 : x^8+1
17 : x^16+x^15+x^14+x^13+x^12+x^11+x^10+x^9+x^8+x^7+x^6+x^5+x^4+x^3+x^2+x+1
18 : x^6-x^3+1
19 : x^18+x^17+x^16+x^15+x^14+x^13+x^12+x^11+x^10+x^9+x^8+x^7+x^6+x^5+x^4+x^3+x^2+x+1
20 : x^8-x^6+x^4-x^2+1

円分多項式の係数は、上の例を見ると、必ず、−１、０、１となるように見えるが、
ｎ＝105 のとき、２となる（105＝３×５×７に注意）。
実は、円分多項式の係数は、任意の整数値を取りうる。