[ < ] [ > ]   [ << ] [ Up ] [ >> ]         [Top] [Contents] [Index] [ ? ]

29. Number Theory


[ < ] [ > ]   [ << ] [ Up ] [ >> ]         [Top] [Contents] [Index] [ ? ]

29.1 Functions and Variables for Number Theory

関数: bern (n)

整数 nについて n番目の Bernoulli数を返します。 もし zerobernfalseならゼロに等しい Bernoulli数は抑制されます。

burnも参照してください。

 
(%i1) zerobern: true$
(%i2) map (bern, [0, 1, 2, 3, 4, 5, 6, 7, 8]);
                      1  1       1      1        1
(%o2)           [1, - -, -, 0, - --, 0, --, 0, - --]
                      2  6       30     42       30
(%i3) zerobern: false$
(%i4) map (bern, [0, 1, 2, 3, 4, 5, 6, 7, 8]);
                      1  1    1   1     1   5     691   7
(%o4)           [1, - -, -, - --, --, - --, --, - ----, -]
                      2  6    30  42    30  66    2730  6

Categories:  Number theory

関数: bernpoly (x, n)

変数 xに関する n番目の Bernoulli多項式を返します。

Categories:  Number theory

関数: bfzeta (s, n)

引数 sに関する Riemannのゼータ関数を返します。 戻り値は多倍長浮動小数点です; nは戻り値の小数点以下の桁数です。

関数: bfhzeta (s, h, n)

引数 shに関する Hurwitzのゼータ関数を返します。 戻り値は多倍長浮動小数点です; nは戻り値の小数点以下の桁数です。

Hurwitzゼータ関数は以下のように定義されます。

 
                        inf
                        ====
                        \        1
         zeta (s,h)  =   >    --------
                        /            s
                        ====  (k + h)
                        k = 0

load ("bffac")でこの関数をロードします。

関数: burn (n)

n番目の Bernoulli数の近似の有理数をを返します。 burnは(有理) Bernoulli数が まあまあの効率で(超越的)ゼータによって近似できるという観察を利用します。

 
                   n - 1  1 - 2 n
              (- 1)      2        zeta(2 n) (2 n)!
     B(2 n) = ------------------------------------
                                2 n
                             %pi

bernは返す前にインデックス nまでの Bernoulli数すべてを計算するので、 burnは大きな、孤立した n(たぶん 105以上の n)に対しては bernより効率的かもしれません。 burnは 255よりおおきな偶数 nに対しては近似を呼び出します。 奇数と 255以下の nに対しては関数 bernが呼び出されます。

load ("bffac")でこの関数をロードします。 bernも参照してください。

Categories:  Number theory

関数: chinese ([r_1, …, r_n], [m_1, …, m_n])

連立合同式 x = r_1 mod m_1, …, x = r_n mod m_nを解きます。 剰余 r_nは任意の整数が可能である一方、 法 m_nは正の互いに素な整数でなければいけません。

 
(%i1) mods : [1000, 1001, 1003, 1007];
(%o1)                   [1000, 1001, 1003, 1007]
(%i2) lreduce('gcd, mods);
(%o2)                               1
(%i3) x : random(apply("*", mods));
(%o3)                         685124877004
(%i4) rems : map(lambda([z], mod(x, z)), mods);
(%o4)                       [4, 568, 54, 624]
(%i5) chinese(rems, mods);
(%o5)                         685124877004
(%i6) chinese([1, 2], [3, n]);
(%o6)                    chinese([1, 2], [3, n])
(%i7) %, n = 4;
(%o7)                              10

Categories:  Number theory

関数: cf (expr)

連分数近似を計算します。 exprは連分数と整数の平方根と実数リテラル(整数、有理数、通常の浮動小数点や多倍長浮動小数点)から構成される式です。 cfは有理数に関して厳密な展開を計算しますが、 通常の浮動小数点に関して展開は ratepsilonで丸められ、 多倍長小数点に関して 10^(-fpprec)で丸められます。

式の中のオペランドは代数演算子を組み合わせられます。 Maximaは cfの外側で連分数に関する演算について知りません。

cfは、 listarithfalseにバインドした後、引数を評価します。 cfはリストとして表現された連分数を返します。

連分数 a + 1/(b + 1/(c + ...))はリスト [a, b, c, ...]で表現されます。 リストの要素 a, b, c, ...は整数に評価されなければいけません。 exprsqrt (n)も含むかもしれません。 nは整数です。 この場合、cfは変数 cflengthの値掛ける周期と同じ数の連分数の項を与えます。

cfdisrepが返す代数表現を評価することで、連分数は数に評価することができます。 連分数を評価する別の方法に関しては cfexpandも参照してください。

cfdisrep, cfexpand, cflengthも参照してください。

例:

Categories:  Continued fractions

関数: cfdisrep (list)

連分数 [a, b, c, ...]のリスト表現から形式 a + 1/(b + 1/(c + ...))の通常の代数式を構成し返します。

 
(%i1) cf ([1, 2, -3] + [1, -2, 1]);
(%o1)                     [1, 1, 1, 2]
(%i2) cfdisrep (%);
                                  1
(%o2)                     1 + ---------
                                    1
                              1 + -----
                                      1
                                  1 + -
                                      2

Categories:  Continued fractions

関数: cfexpand (x)

連分数 xのコンバージェントの最後(列1)とその1つ前(列2)の分子と分母の行列を返します。

 
(%i1) cf (rat (ev (%pi, numer)));

`rat' replaced 3.141592653589793 by 103993/33102 =3.141592653011902
(%o1)                  [3, 7, 15, 1, 292]
(%i2) cfexpand (%);
                         [ 103993  355 ]
(%o2)                    [             ]
                         [ 33102   113 ]
(%i3) %[1,1]/%[2,1], numer;
(%o3)                   3.141592653011902

Categories:  Continued fractions

オプション変数: cflength

デフォルト値: 1

cflengthは、値 cflength掛ける周期として関数 cfが与える連分数の項の数を制御します。 従って、デフォルトは 1周期を与えます。

 
(%i1) cflength: 1$
(%i2) cf ((1 + sqrt(5))/2);
(%o2)                    [1, 1, 1, 1, 2]
(%i3) cflength: 2$
(%i4) cf ((1 + sqrt(5))/2);
(%o4)               [1, 1, 1, 1, 1, 1, 1, 2]
(%i5) cflength: 3$
(%i6) cf ((1 + sqrt(5))/2);
(%o6)           [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2]

Categories:  Continued fractions

関数: divsum  
    divsum (n, k)  
    divsum (n)

divsum (n, k)nの約数の k乗した和を返します。

divsum (n)nの約数の和を返します。

 
(%i1) divsum (12);
(%o1)                          28
(%i2) 1 + 2 + 3 + 4 + 6 + 12;
(%o2)                          28
(%i3) divsum (12, 2);
(%o3)                          210
(%i4) 1^2 + 2^2 + 3^2 + 4^2 + 6^2 + 12^2;
(%o4)                          210

Categories:  Number theory

関数: euler (n)

非負の整数 nに対して n番目のEuler数を返します。 もし zerobernfalseなら、 0のEuler数は抑制されます。

Euler-Mascheroni定数に関しては %gammaを参照してください。

 
(%i1) zerobern: true$
(%i2) map (euler, [0, 1, 2, 3, 4, 5, 6]);
(%o2)               [1, 0, - 1, 0, 5, 0, - 61]
(%i3) zerobern: false$
(%i4) map (euler, [0, 1, 2, 3, 4, 5, 6]);
(%o4)               [1, - 1, 5, - 61, 1385, - 50521, 2702765]

Categories:  Number theory

オプション変数: factors_only

デフォルト値: false

ifactorsが返す値を制御します。 デフォルトの falseの時、 ifactorsは計算された素因数の多重性について情報を提供するようになります。 もし factors_onlytrueに設定されているなら、 ifactorsは素因数のリスト以外なにも返しません。

例: ifactorsを参照してください。

Categories:  Number theory

関数: fib (n)

n項の Fibonacci数を返します。 fib(0)は 0に等しく、 fib(1)は 1に等しく、 fib (-n)(-1)^(n + 1) * fib(n)に等しいです。

fibをコールした後, prevfibは 最後に計算された1つ前の Fibonacci数 fib (x - 1)に等しいです。

 
(%i1) map (fib, [-4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8]);
(%o1)           [- 3, 2, - 1, 1, 0, 1, 1, 2, 3, 5, 8, 13, 21]

Categories:  Number theory

関数: fibtophi (expr)

exprに関する Fibonacci数を定数 %phiを使って表現します。 %phi(1 + sqrt(5))/2, 近似的に1.61803399です。

例:

 
(%i1) fibtophi (fib (n));
                           n             n
                       %phi  - (1 - %phi)
(%o1)                  -------------------
                           2 %phi - 1
(%i2) fib (n-1) + fib (n) - fib (n+1);
(%o2)          - fib(n + 1) + fib(n) + fib(n - 1)
(%i3) fibtophi (%);
            n + 1             n + 1       n             n
        %phi      - (1 - %phi)        %phi  - (1 - %phi)
(%o3) - --------------------------- + -------------------
                2 %phi - 1                2 %phi - 1
                                          n - 1             n - 1
                                      %phi      - (1 - %phi)
                                    + ---------------------------
                                              2 %phi - 1
(%i4) ratsimp (%);
(%o4)                           0

Categories:  Number theory

関数: ifactors (n)

正の整数 nに対して nの素因数分解を返します。 もし n=p1^e1..pk^nknの素因数への分解なら、 ifactorsは [[p1, e1], ... , [pk, ek]]を返します。

使われる素因数分解法は 9973までの素数による試行除算と、 Pollardのローとp-1法と、楕円曲線法です。

もし変数 ifactor_verbosetrueに設定されているなら、 ifactorは、因子が見つかり次第直近のフィードバックを含めて何をしているか詳細出力を生成します。

ifactorsが返す値はオプション変数 factors_onlyによって制御されます。 デフォルトの falseの時、 ifactorsは計算された素因数の多重性について情報を提供するようになります。 もし factors_onlytrueに設定されているなら、 ifactorsは単に素因数のリストを返します。

 
(%i1) ifactors(51575319651600);
(%o1)     [[2, 4], [3, 2], [5, 2], [1583, 1], [9050207, 1]]
(%i2) apply("*", map(lambda([u], u[1]^u[2]), %));
(%o2)                        51575319651600
(%i3) ifactors(51575319651600), factors_only : true;
(%o3)                   [2, 3, 5, 1583, 9050207]

Categories:  Number theory

関数: igcdex (n, k)

リスト [a, b, u]を返します。 ここで、 unkの最大公約数で、 ua n + b kに等しいです。 引数 nkは整数でなければいけません。

igcdexはユークリッドのアルゴリズムを実装します。 gcdex.も参照してください。

コマンド load(gcdex)でこの関数をロードします。

例:

 
(%i1) load(gcdex)$

(%i2) igcdex(30,18);
(%o2)                      [- 1, 2, 6]
(%i3) igcdex(1526757668, 7835626735736);
(%o3)            [845922341123, - 164826435, 4]
(%i4) igcdex(fib(20), fib(21));
(%o4)                   [4181, - 2584, 1]

Categories:  Number theory

関数: inrt (x, n)

xの絶対値の整数 n乗根を返します。

 
(%i1) l: [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]$
(%i2) map (lambda ([a], inrt (10^a, 3)), l);
(%o2) [2, 4, 10, 21, 46, 100, 215, 464, 1000, 2154, 4641, 10000]

Categories:  Number theory

関数: inv_mod (n, m)

mを法とする nの逆元を計算します。 もし nmを法とするゼロ因子なら、 inv_mod (n,m)falseを返します。

 
(%i1) inv_mod(3, 41);
(%o1)                           14
(%i2) ratsimp(3^-1), modulus = 41;
(%o2)                           14
(%i3) inv_mod(3, 42);
(%o3)                          false

Categories:  Number theory

関数: isqrt (x)

整数 xの絶対値の「整数平方根」を返します。

Categories:  Mathematical functions

関数: jacobi (p, q)

pqの Jacobi記号を返します。

 
(%i1) l: [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]$
(%i2) map (lambda ([a], jacobi (a, 9)), l);
(%o2)         [1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0]

Categories:  Number theory

関数: lcm (expr_1, ..., expr_n)

引数の最小公倍数を返します。 引数は、整数はもちろん一般式を取り得ます。

load ("functs")でこの関数をロードします。

Categories:  Number theory

関数: lucas (n)

n番目の Lucas数を返します。 lucas(0)は 2に等しく、 lucas(1)は 1に等しく、 lucas(-n)(-1)^(-n) * lucas(n)に等しいです。

 
(%i1) map (lucas, [-4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8]);
(%o1)             [7, - 4, 3, - 1, 2, 1, 3, 4, 7, 11, 18, 29, 47]

After calling lucasを呼び出した後には、グローバル変数 next_lucasは 最後の戻り値に続く Lucas数 lucas (n + 1)に等しいです。 例は、どのようにして Fibonacci数が lucasnext_lucasを使って計算可能かを示します。

 
(%i1) fib_via_lucas(n) :=
         block([lucas : lucas(n)],
         signum(n) * (2*next_lucas - lucas)/5 )$
(%i2) map (fib_via_lucas, [-4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8]);
(%o2)             [- 3, 2, - 1, 1, 0, 1, 1, 2, 3, 5, 8, 13, 21]

Categories:  Number theory

関数: mod (x, y)

もし xyが実数で、 yがゼロでないなら、 x - y * floor(x / y)を返します。 さらにすべての実数 xに関して mod (x, 0) = xが成り立ちます。 定義 mod (x, 0) = xの議論に関しては、 Graham, Knuth, Patashnik著の「コンピュータの数学」の3.4節を参照してください。 関数 mod (x, 1) は周期が 1で mod (1, 1) = 0mod (0, 1) = 0ののこぎり波関数です。

複素数の偏角の主値(区間 (-%pi, %pi]での数)を見つけるためには、 関数 x |-> %pi - mod (%pi - x, 2*%pi)を使います。 xは引数です。

xyが定数式(例えば 10 * %pi)の時、 modfloorceilingが使うのと同じ多倍長浮動小数点評価スキームを使います。 同様に、まれですが そんな場合にはmodは間違った値を返すことがありえます。

数値でない引数 xyに関して modはいくつかの式整理規則を知っています:

 
(%i1) mod (x, 0);
(%o1)                           x
(%i2) mod (a*x, a*y);
(%o2)                      a mod(x, y)
(%i3) mod (0, x);
(%o3)                           0

Categories:  Mathematical functions

関数: next_prime (n)

nよりも大きな最も小さな素数を返します。

 
(%i1) next_prime(27);
(%o1)                       29

Categories:  Number theory

関数: partfrac (expr, var)

主変数 varに関する部分分数式 exprを展開します。 partfracは完全な部分分数分解を行います。 利用したアルゴリズムは 部分分数展開(元の分母の因子)の分母は互いに素であるという事実に基づいています。 分子は分母の線形結合として書けて、結果は展開されたものになります。

 
(%i1) 1/(1+x)^2 - 2/(1+x) + 2/(2+x);
                      2       2        1
(%o1)               ----- - ----- + --------
                    x + 2   x + 1          2
                                    (x + 1)
(%i2) ratsimp (%);
                                 x
(%o2)                 - -------------------
                         3      2
                        x  + 4 x  + 5 x + 2
(%i3) partfrac (%, x);
                      2       2        1
(%o3)               ----- - ----- + --------
                    x + 2   x + 1          2
                                    (x + 1)

関数: power_mod (a, n, m)

a^n mod mを計算するために剰余アルゴリズムを使います。 ここで、 anは整数で、 mは正の整数です。 もし nが負なら inv_modが剰余逆元を見つけるために使われます。

 
(%i1) power_mod(3, 15, 5);
(%o1)                          2
<(%i2) mod(3^15,5);
(%o2)                          2
(%i3) power_mod(2, -1, 5);
(%o3)                          3
(%i4) inv_mod(2,5);
(%o4)                          3

Categories:  Number theory

関数: primep (n)

素数テスト。 もし primep (n)falseを返すなら、 nは合成数であり、 もし trueを返すなら、 nは非常に高い確率で素数です。

341550071728321より小さな nに対しては Miller-Rabinのテストの決定的バージョンが使われます。 もし primep (n)trueを返すなら、 nは素数です。

341550071728321よりの大きな nに対して、 primepは、 primep_number_of_tests個の Miller-Rabinの疑似素数テストと 1つの Lucasの疑似素数テストを使います。 合成数の nが Miller-Rabinのテスト1つを通過する確率は 1/4より小さいです。 primep_number_of_testsに関してデフォルト値 25を使うと、 通過した nが合成である確率は 10^-15よりもはるかに小さいです。

オプション変数: primep_number_of_tests

デフォルト値: 25

primepの中で使われる Miller-Rabinのテストの回数。

Categories:  Number theory

関数: primes (start, end)

startから endまでのすべての素数のリストを返します。

 
(%i1) primes(3, 7);
(%o1)                     [3, 5, 7]

Categories:  Number theory

関数: prev_prime (n)

nよりも小さな最大の素数を返します。

 
(%i1) prev_prime(27);
(%o1)                       23

Categories:  Number theory

関数: qunit (n)

実二次数体 sqrt (n)の基本単数、すなわちノルムが 1の要素を返します。 ここで nは整数です。 これは結果的にペル方程式 a^2 - n b^2 = 1を解くことになります。

 
(%i1) qunit (17);
(%o1)                     sqrt(17) + 4
(%i2) expand (% * (sqrt(17) - 4));
(%o2)                           1

Categories:  Number theory

関数: totient (n)

n以下の nと互いに素な整数の数を返します。

Categories:  Number theory

オプション変数: zerobern

デフォルト値: true

zerobernfalseの時、 bernは Bernoulli数を除外し、 eulerはゼロに等しい Euler数を除外します。 berneulerを参照してください。

Categories:  Number theory

関数: zeta (n)

Riemannのゼータ関数を返します。 もし xが負の整数か, 0, 1,または正の偶数なら、 Reimannのゼータ関数は厳密な値に整理されます。 正の偶数に対してはオプション変数 zeta%pitrueであることも必要です。 (zeta%piを参照してください。) 浮動小数点または多倍長浮動小数点数に対しては Reimannゼータ関数は数値的に評価されます。 Maximaは、 有理非整数、浮動小数点数、複素数の引数を含む他の引数すべてに対して、 また、 zeta%piが値 falseなら偶数に対しても 名詞形 zeta (n)を返します。

zeta(1)は未定義ですが、 Maximaは上からと下からの極限 limit(zeta(x), x, ,1)を知っています。

bfzetazeta%piも参照してください。

例:

 
(%i1) zeta([-2, -1, 0, 0.5, 2, 3, 1+%i]);
                                             2
            1     1                       %pi
(%o1) [0, - --, - -, - 1.460354508809586, ----, zeta(3),
            12    2                        6
                                                    zeta(%i + 1)]
(%i2) limit(zeta(x),x,1,plus);
(%o2)                          inf
(%i3) limit(zeta(x),x,1,minus);
(%o3)                         minf

Categories:  Number theory

オプション変数: zeta%pi

デフォルト値: true

zeta%pitrueの時、 偶数 nに対して zeta%pi^nに比例する式を返します。 そうでないなら、 偶数 nに対して zetaは名詞形 zeta (n)を返します。

例:

 
(%i1) zeta%pi: true$
(%i2) zeta (4);
                                 4
                              %pi
(%o2)                         ----
                               90
(%i3) zeta%pi: false$
(%i4) zeta (4);
(%o4)                        zeta(4)

Categories:  Number theory

関数: zn_add_table (n)

(Z/nZ)のすべての要素の加算表を表示します。

zn_mult_table, zn_power_tableも参照してください。

Categories:  Number theory

関数: zn_characteristic_factors (n)

nのトーティエントの特性因子を含むリストを返します。

特性因子を使って、 nを法とする乗法群を巡回部分群の群直積として表現できます。

群自身が巡回的の時には、リストはトーティエントのみを含み、 zn_primrootを使って使って生成元を計算できます。 もしトーティエントが複数の特性因子に分割されるなら、 zn_factor_generatorsは対応する部分群の生成元を見つけます。

リストの r個の因子のそれぞれは右に続く因子を割り切ります。 従って、最後の因子 f_rnと互いに素なすべての aに対して a^f_r = 1 (mod n)を満たします。 この因子は Carmichael函数や Carmichaelラムダとしても知られています。

もし n > 2なら、 totient(n)/2^rは平方剰余の数であり、 これらのそれぞれは 2^r個の平方根を持ちます。

totient, zn_primroot, zn_factor_generatorsも参照してください。

例:

14を法とする乗法群は巡回的で、その 6要素は原始根で生成できます。

 
(%i1) [zn_characteristic_factors(14), phi: totient(14)];
(%o1)                              [[6], 6]
(%i2) [zn_factor_generators(14), g: zn_primroot(14)];
(%o2)                              [[3], 3]
(%i3) M14: makelist(power_mod(g,i,14), i,0,phi-1);
(%o3)                         [1, 3, 9, 13, 11, 5]

15を法とする乗法群は巡回的でなかく、その 8要素は2つの因子生成元で生成できます。

 
(%i1) [[f1,f2]: zn_characteristic_factors(15), totient(15)];
(%o1)                             [[2, 4], 8]
(%i2) [[g1,g2]: zn_factor_generators(15), zn_primroot(15)];
(%o2)                           [[11, 7], false]
(%i3) UG1: makelist(power_mod(g1,i,15), i,0,f1-1);
(%o3)                               [1, 11]
(%i4) UG2: makelist(power_mod(g2,i,15), i,0,f2-1);
(%o4)                            [1, 7, 4, 13]
(%i5) M15: create_list(mod(i*j,15), i,UG1, j,UG2);
(%o5)                      [1, 7, 4, 13, 11, 2, 14, 8]

最後の特性因子 4に関して、 M15の中のすべての aに対して a^4 = 1 (mod 15)を満たします。

M15は2つの特性因子と、8/2^2平方剰余を持ち、 これらのそれぞれは 2^2個の平方根を持ちます。

 
(%i6) zn_power_table(15);
                               [ 1   1  1   1 ]
                               [              ]
                               [ 2   4  8   1 ]
                               [              ]
                               [ 4   1  4   1 ]
                               [              ]
                               [ 7   4  13  1 ]
(%o6)                          [              ]
                               [ 8   4  2   1 ]
                               [              ]
                               [ 11  1  11  1 ]
                               [              ]
                               [ 13  4  7   1 ]
                               [              ]
                               [ 14  1  14  1 ]
(%i7) map(lambda([i], zn_nth_root(i,2,15)), [1,4]);
(%o7)                   [[1, 4, 11, 14], [2, 7, 8, 13]]

Categories:  Number theory

関数: zn_carmichael_lambda (n)

もし n1なら 1を返し、 そうでなければ、 nのトーシェントの最大特性因子を返します。

注意と例は zn_characteristic_factorsを参照してください。

Categories:  Number theory

関数: zn_determinant (matrix, p)

LU分解の技法を使って、 (Z/pZ)上の matrixの行列式を計算します。 pは素数でなければいけません。

行列式がゼロに等しい場合、LU分解が失敗するかもしれません。 この場合、 zn_determinantはモジュラーでない行列式を計算し、その後整理します。

zn_invert_by_luも参照してください。

例:

 
(%i1) m : matrix([1,3],[2,4]);
                                [ 1  3 ]
(%o1)                           [      ]
                                [ 2  4 ]
(%i2) zn_determinant(m, 5);
(%o2)                               3
(%i3) m : matrix([2,4,1],[3,1,4],[4,3,2]);
                               [ 2  4  1 ]
                               [         ]
(%o3)                          [ 3  1  4 ]
                               [         ]
                               [ 4  3  2 ]
(%i4) zn_determinant(m, 5);
(%o4)                               0

Categories:  Number theory

関数: zn_factor_generators (n)

nのトーティエントの特性因子に対応する因子生成元を含むリストを返します。

コメントと例に関しては zn_characteristic_factorsを参照してください。

Categories:  Number theory

関数: zn_invert_by_lu (matrix, p)

LU分解の技法を使って、 (Z/pZ)上で matrixのモジュラー逆元を計算します。 pは素数、 matrixは可逆でなければいけません。 matrixが可逆でないなら、 zn_invert_by_lufalseを返します。

zn_determinantも参照してください。

例:

 
(%i1) m : matrix([1,3],[2,4]);
                                [ 1  3 ]
(%o1)                           [      ]
                                [ 2  4 ]
(%i2) zn_determinant(m, 5);
(%o2)                               3
(%i3) mi : zn_invert_by_lu(m, 5);
                                [ 3  4 ]
(%o3)                           [      ]
                                [ 1  2 ]
(%i4) matrixmap(lambda([a], mod(a, 5)), m . mi);
                                [ 1  0 ]
(%o4)                           [      ]
                                [ 0  1 ]

Categories:  Number theory

関数: zn_log  
    zn_log (a, g, n)  
    zn_log (a, g, n, [[p1, e1], …, [pk, ek]])

離散対数を計算します。 (Z/nZ)* が巡回群、 gnを法とする原始根とし、 aがこの群の要素とします。 そのとき、zn_log (a, g, n)は合同式 g^x = a mod nを解きます。

採用したアルゴリズムは totient(n)の素因数分解を必要とします。 この素因数分解はその上時間を消費するかもしれず、いくつかの場合、 最初に素因数分解してそれから因数のリストを 4番目の引数として zn_logに渡すのが実用的かもしれません。 リストは、デフォルトオプション factors_only : falseを使って ifactors(totient(n))が返すリストと同じ形式でなければいけません。

離散対数のために、アルゴリズムは Pohlig-Hellman-縮約と Pollardのロー法を使います。 zn_logの実行時間は主にtotientの最大素因数のビット長に依存します。

zn_primroot, zn_order, ifactors, totientも参照して下さい。

例:

zn_log (a, g, n)は合同式 g^x = a mod nを解きます。

 
(%i1) n : 22$
(%i2) g : zn_primroot(n);
(%o2)                               7
(%i3) ord_7 : zn_order(7, n);
(%o3)                              10
(%i4) powers_7 : makelist(power_mod(g, x, n), x, 0, ord_7 - 1);
(%o4)              [1, 7, 5, 13, 3, 21, 15, 17, 9, 19]
(%i5) zn_log(21, g, n);
(%o5)                               5
(%i6) map(lambda([x], zn_log(x, g, n)), powers_7);
(%o6)                [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]

オプションの4番目の引数は ifactors(totient(n))が返すリストと同じ形式でなければいけません。 実行時間は主にtotientの最大素因数のビット長に依存します。

 
(%i1) (p : 2^127-1, primep(p));
(%o1)                             true
(%i2) ifs : ifactors(p - 1)$
(%i3) g : zn_primroot(p, ifs);
(%o3)                              43
(%i4) a : power_mod(g, 1234567890, p)$
(%i5) zn_log(a, g, p, ifs);
(%o5)                          1234567890
(%i6) time(%o5);
(%o6)                            [1.204]
(%i7) f_max : last(ifs);
(%o7)                       [77158673929, 1]
(%i8) slength( printf(false, "~b", f_max[1]) );
(%o8)                              37

Categories:  Number theory

関数: zn_mult_table  
    zn_mult_table (n)  
    zn_mult_table (n, gcd)

オプション引数 gcdがない時、 zn_mult_table(n)は、 nと互いに素な (Z/nZ)*のすべての要素の乗算表を表示します。

オプションの第二引数 gcdは、(Z/nZ)の特定のぶぶ集合を選択することを可能にします。 もし gcdが整数なら、 gcd(x,n) = gcdであるすべての剰余 xの乗算表を返します。 読みやすくするため行と列の見出しも追加されます。 必要なら、これらは submatrix(1, table, 1)で簡単に取り除くことができます。

もし gcdとして allを渡したら、 (Z/nZ)のゼロでないすべての要素を表に印字します。

2番目の例は部分群の乗算表を生成する代わりの方法を提示します。

zn_add_table, zn_power_tableも参照してください。

例:

デフォルトの表は (Z/nZ)*のすべての要素を表示し、 これでモジュラー乗法群の基本性質を例示したり勉強したりできます。 例えば、主対角線はすべての平方剰余を含み、 それぞれの行と列はすべての要素を含み、表は対象であるなど。

もし gcdallを与えると、(Z/nZ)のゼロでない要素すべてを表で表示します。

 
(%i1) zn_mult_table(8);
                                [ 1  3  5  7 ]
                                [            ]
                                [ 3  1  7  5 ]
(%o1)                           [            ]
                                [ 5  7  1  3 ]
                                [            ]
                                [ 7  5  3  1 ]
(%i2) zn_mult_table(8, all);
                            [ 1  2  3  4  5  6  7 ]
                            [                     ]
                            [ 2  4  6  0  2  4  6 ]
                            [                     ]
                            [ 3  6  1  4  7  2  5 ]
                            [                     ]
(%o2)                       [ 4  0  4  0  4  0  4 ]
                            [                     ]
                            [ 5  2  7  4  1  6  3 ]
                            [                     ]
                            [ 6  4  2  0  6  4  2 ]
                            [                     ]
                            [ 7  6  5  4  3  2  1 ]

もし gcdが整数なら、 見やすくするため行と列の見出しがつきます。

もし gcdで選ばれた部分集合が群なら、乗算表を生成する別の方法があります。 1を単位元として持つ群からの同型写像は読みやすい表を作ります。 写像はCRTを介して行われます。

T36_4の第2版の中で、ちょうど表 T9のように、単位元、ここでは 28、が左上隅に置かれます。

 
(%i1) T36_4: zn_mult_table(36,4);
                        [ *   4   8   16  20  28  32 ]
                        [                            ]
                        [ 4   16  32  28  8   4   20 ]
                        [                            ]
                        [ 8   32  28  20  16  8   4  ]
                        [                            ]
(%o1)                   [ 16  28  20  4   32  16  8  ]
                        [                            ]
                        [ 20  8   16  32  4   20  28 ]
                        [                            ]
                        [ 28  4   8   16  20  28  32 ]
                        [                            ]
                        [ 32  20  4   8   28  32  16 ]
(%i2) T9: zn_mult_table(36/4);
                             [ 1  2  4  5  7  8 ]
                             [                  ]
                             [ 2  4  8  1  5  7 ]
                             [                  ]
                             [ 4  8  7  2  1  5 ]
(%o2)                        [                  ]
                             [ 5  1  2  7  8  4 ]
                             [                  ]
                             [ 7  5  1  8  4  2 ]
                             [                  ]
                             [ 8  7  5  4  2  1 ]
(%i3) T36_4: matrixmap(lambda([x], chinese([0,x],[4,9])), T9);
                          [ 28  20  4   32  16  8  ]
                          [                        ]
                          [ 20  4   8   28  32  16 ]
                          [                        ]
                          [ 4   8   16  20  28  32 ]
(%o3)                     [                        ]
                          [ 32  28  20  16  8   4  ]
                          [                        ]
                          [ 16  32  28  8   4   20 ]
                          [                        ]
                          [ 8   16  32  4   20  28 ]

Categories:  Number theory

関数: zn_nth_root  
    zn_nth_root (x, n, m)  
    zn_nth_root (x, n, m, [[p1, e1], …, [pk, ek]])

xを含む (Z/mZ)の乗法部分群から xn乗根すべてのリストを返します。 もし xが法 mに対する n冪剰余でないか、 (Z/mZ)のどの乗法部分群にも含まれないなら、 falseを返します。

もし最大公約数 g = gcd(x,m)m/gと互いに素なら、 xは法 mに対する乗法部分群の要素です。

zn_nth_rootは Adlemanと Manders、 Millerによるアルゴリズムと Daniel Shanksによるモジュラー乗法群に関する定理に基づいています。

アルゴリズムはモジュラス mの素因数分解を必要とします。 mの因数分解がわかっている場合、4番目の引数として因子のリストを渡すことができます。 このオプション引数は デフォルトオプション factors_only: falseを使った時 ifactors(m)が返すリストと同じ形式のものでなければいけません。

例:

n1から6までの値をとる、 1n番目の根すべてを含むリストのリストが続く 14を法とする乗法群のべきテーブル

 
(%i1) zn_power_table(14);
                         [ 1   1   1   1   1   1 ]
                         [                       ]
                         [ 3   9   13  11  5   1 ]
                         [                       ]
                         [ 5   11  13  9   3   1 ]
(%o1)                    [                       ]
                         [ 9   11  1   9   11  1 ]
                         [                       ]
                         [ 11  9   1   11  9   1 ]
                         [                       ]
                         [ 13  1   13  1   13  1 ]
(%i2) makelist(zn_nth_root(1,n,14), n,1,6);
(%o2)  [[1], [1, 13], [1, 9, 11], [1, 13], [1], [1, 3, 5, 9, 11, 13]]

以下の例では xmと互いに素ではないですが、 (Z/mZ)の乗法部分群の要素であり、 任意の n乗根は同じ部分群の要素です。

剰余クラス 3は どの (Z/63Z)の乗法部分群の要素でもなく、 なので 27の3乗根として返されません。

ここで zn_power_tablegcd(x,63) = 9が成り立つ (Z/63Z)のすべての剰余 xを表示します。 この部分群は (Z/7Z)*に同型で、単位元 36はCRTを介して計算されます。

 
(%i1) m: 7*9$

(%i2) zn_power_table(m,9);
                         [ 9   18  36  9   18  36 ]
                         [                        ]
                         [ 18  9   36  18  9   36 ]
                         [                        ]
                         [ 27  36  27  36  27  36 ]
(%o2)                    [                        ]
                         [ 36  36  36  36  36  36 ]
                         [                        ]
                         [ 45  9   27  18  54  36 ]
                         [                        ]
                         [ 54  18  27  9   45  36 ]
(%i3) zn_nth_root(27,3,m);
(%o3)                           [27, 45, 54]
(%i4) id7:1$  id63_9: chinese([id7,0],[7,9]);
(%o5)                                36

以下のRSAのような例では、モジュラス Nは平方因子をもちません。 すなわち、排他的に1乗冪因子に分解され、 0から N-1までのすべての xは乗法部分群に含まれます。

解読処理は e乗根を必要とします。 etotient(N)と互いに素で、なので e乗根は一意です。 この場合、 zn_nth_rootは効果的に CRT-RSAを実行します。 (flattenは括弧は取り除きますが解は取り除かないことに注意してください)

 
(%i1) [p,q,e]: [5,7,17]$  N: p*q$

(%i3) xs: makelist(x,x,0,N-1)$

(%i4) ys: map(lambda([x],power_mod(x,e,N)),xs)$

(%i5) zs: flatten(map(lambda([y], zn_nth_root(y,e,N)), ys))$

(%i6) is(zs = xs);
(%o6)                             true

以下の例では、モジュラスの因数分解が知られていて、 4番目の引数として渡されます。

 
(%i1) p: 2^107-1$  q: 2^127-1$  N: p*q$

(%i4) ibase: obase: 16$

(%i5) msg: 11223344556677889900aabbccddeeff$

(%i6) enc: power_mod(msg, 10001, N);
(%o6)    1a8db7892ae588bdc2be25dd5107a425001fe9c82161abc673241c8b383
(%i7) zn_nth_root(enc, 10001, N, [[p,1],[q,1]]);
(%o7)               [11223344556677889900aabbccddeeff]

Categories:  Number theory

関数: zn_order  
    zn_order (x, n)  
    zn_order (x, n, [[p1, e1], …, [pk, ek]])

xが 有限群 (Z/nZ)*の単位元ならその次数を返し、そうでないなら falseを返します。 xnと互いに素なら nを法とする単位元です。

採用したアルゴリズムは totient(n)の素因数分解を必要とします。 この素因数分解はその上時間を消費するかもしれず、いくつかの場合、 最初に素因数分解してそれから因数のリストを 3番目の引数として zn_orderに渡すのが実用的かもしれません。 リストは、デフォルトオプション factors_only : falseを使って ifactors(totient(n))が返すリストと同じ形式でなければいけません。

zn_primroot, ifactors, totientも参照して下さい。

例:

zn_orderは (Z/nZ)*の単位元 xの次数を計算します。

 
(%i1) n : 22$
(%i2) g : zn_primroot(n);
(%o2)                               7
(%i3) units_22 : sublist(makelist(i,i,1,21), lambda([x], gcd(x, n) = 1));
(%o3)              [1, 3, 5, 7, 9, 13, 15, 17, 19, 21]
(%i4) (ord_7 : zn_order(7, n)) = totient(n);
(%o4)                            10 = 10
(%i5) powers_7 : makelist(power_mod(g,i,n), i,0,ord_7 - 1);
(%o5)              [1, 7, 5, 13, 3, 21, 15, 17, 9, 19]
(%i6) map(lambda([x], zn_order(x, n)), powers_7);
(%o6)              [1, 10, 5, 10, 5, 2, 5, 10, 5, 10]
(%i7) map(lambda([x], ord_7/gcd(x, ord_7)), makelist(i, i,0,ord_7 - 1));
(%o7)              [1, 10, 5, 10, 5, 2, 5, 10, 5, 10]
(%i8) totient(totient(n));
(%o8)                               4

オプションの三番目の引数は ifactors(totient(n))が返すリストと同じ形式でなければいけません。

 
(%i1) (p : 2^142 + 217, primep(p));
(%o1)                             true
(%i2) ifs : ifactors( totient(p) )$
(%i3) g : zn_primroot(p, ifs);
(%o3)                               3
(%i4) is( (ord_3 : zn_order(g, p, ifs)) = totient(p) );
(%o4)                             true
(%i5) map(lambda([x], ord_3/zn_order(x, p, ifs)), makelist(i,i,2,15));
(%o5)        [22, 1, 44, 10, 5, 2, 22, 2, 8, 2, 1, 1, 20, 1]

Categories:  Number theory

関数: zn_power_table  
    zn_power_table (n)  
    zn_power_table (n, gcd)  
    zn_power_table (n, gcd, max_exp)

オプション引数なしの場合、 zn_power_table(n)nと互いに素なすべての剰余クラスである (Z/nZ)*のすべての要素のべき表を表示します。 指数は 1から(Carmichael函数や Carmichaelラムダとしても知られる) totient(n) の最大特性因子までループし、表は右側1の列で終わります。

オプションの第2引数 gcdで (Z/nZ)の特定の部分集合の冪を選択することができます。 もし gcdが整数なら、 gcd(x,n) = gcdが成立するすべての剰余クラス xの冪を返します。 すなわち、 gcdのデフォルト値は 1です。 もし gcdallを渡したら、 表は (Z/nZ)のすべての要素の冪を含みます。

もしオプションの第3引数 max_expを与えると、 指数は1から max_expまでループします。

zn_add_table, zn_mult_tableも参照してください。

例:

gcd = 1であるデフォルトでは 例えば FermatやEulerの基本定理を例示したり勉強したりできます。

引数 gcdで (Z/nZ)の部分集合を選択して乗法部分群や同型写像が勉強できます。 例えば、群 G10G10_2は乗算の下でどちらも G5と同型です。 1G5では単位元です。 それぞれで 1そうなります。 G10では 6も単位元です、G10_2 原始根、n乗根などの対応する写像が存在します。

 
(%i1) zn_power_table(10);
                              [ 1  1  1  1 ]
                              [            ]
                              [ 3  9  7  1 ]
(%o1)                         [            ]
                              [ 7  9  3  1 ]
                              [            ]
                              [ 9  1  9  1 ]
(%i2) zn_power_table(10,2);
                              [ 2  4  8  6 ]
                              [            ]
                              [ 4  6  4  6 ]
(%o2)                         [            ]
                              [ 6  6  6  6 ]
                              [            ]
                              [ 8  4  2  6 ]
(%i3) zn_power_table(10,5);
(%o3)                         [ 5  5  5  5 ]
(%i4) zn_power_table(10,10);
(%o4)                         [ 0  0  0  0 ]
(%i5) G5: [1,2,3,4];
(%o6)                          [1, 2, 3, 4]
(%i6) G10_2: map(lambda([x], chinese([0,x],[2,5])), G5);
(%o6)                          [6, 2, 8, 4]
(%i7) G10: map(lambda([x], power_mod(3, zn_log(x,2,5), 10)), G5);
(%o7)                          [1, 3, 7, 9]

もし gcdallを渡せば、 表は (Z/nZ)のすべての要素の冪を含みます。

3番目の引数 max_expで最高指数を設定できます。 以下の表は RSAの非常に小さな例を表示します。

 
(%i1) N:2*5$ phi:totient(N)$ e:7$ d:inv_mod(e,phi)$

(%i5) zn_power_table(N, all, e*d);
       [ 0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 ]
       [                                                               ]
       [ 1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1 ]
       [                                                               ]
       [ 2  4  8  6  2  4  8  6  2  4  8  6  2  4  8  6  2  4  8  6  2 ]
       [                                                               ]
       [ 3  9  7  1  3  9  7  1  3  9  7  1  3  9  7  1  3  9  7  1  3 ]
       [                                                               ]
       [ 4  6  4  6  4  6  4  6  4  6  4  6  4  6  4  6  4  6  4  6  4 ]
(%o5)  [                                                               ]
       [ 5  5  5  5  5  5  5  5  5  5  5  5  5  5  5  5  5  5  5  5  5 ]
       [                                                               ]
       [ 6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6 ]
       [                                                               ]
       [ 7  9  3  1  7  9  3  1  7  9  3  1  7  9  3  1  7  9  3  1  7 ]
       [                                                               ]
       [ 8  4  2  6  8  4  2  6  8  4  2  6  8  4  2  6  8  4  2  6  8 ]
       [                                                               ]
       [ 9  1  9  1  9  1  9  1  9  1  9  1  9  1  9  1  9  1  9  1  9 ]

Categories:  Number theory

関数: zn_primroot  
    zn_primroot (n)  
    zn_primroot (n, [[p1, e1], …, [pk, ek]])

もし乗法群 (Z/nZ)*が巡回的なら、 zn_primrootnを法とする最小の原始根を計算します。 もし n2, 4, p^k または 2*p^kに等しいなら (Z/nZ)* は巡回的です。 ここで pは素数で 2より大きく kは自然数です。 もしオプション変数 zn_primroot_pretest (デフォルト: false)が trueに設定されているなら zn_primrootは条件付き事前テスト(according pretest)を実行します。 どんな場合でも計算は上限 zn_primroot_limitで制限されます。

もし (Z/nZ)*が巡回的でないか zn_primroot_limitまで原始根がないなら、 zn_primrootfalseを返します。

採用したアルゴリズムは totient(n)の素因数分解を必要とします。 この素因数分解はその上時間を消費するかもしれず、いくつかの場合、 最初に素因数分解してそれから因数のリストを 3番目の引数として zn_primrootに渡すのが実用的かもしれません。

リストは、デフォルトオプション factors_only : falseを使って ifactors(totient(n))が返すリストと同じ形式でなければいけません。

zn_primroot_p, zn_order, ifactors, totientも参照してください。

例:

zn_primrootnを法とする最小原始根を計算するか、 falseを返します。

 
(%i1) n : 14$
(%i2) g : zn_primroot(n);
(%o2)                               3
(%i3) zn_order(g, n) = totient(n);
(%o3)                             6 = 6
(%i4) n : 15$
(%i5) zn_primroot(n);
(%o5)                             false

オプションの二番目の引数は ifactors(totient(n))が返すリストと同じ形式でなければいけません。

 
(%i1) (p : 2^142 + 217, primep(p));
(%o1)                             true
(%i2) ifs : ifactors( totient(p) )$
(%i3) g : zn_primroot(p, ifs);
(%o3)                               3
(%i4) [time(%o2), time(%o3)];
(%o4)                    [[15.556972], [0.004]]
(%i5) is(zn_order(g, p, ifs) = p - 1);
(%o5)                             true
(%i6) n : 2^142 + 216$
(%i7) ifs : ifactors(totient(n))$
(%i8) zn_primroot(n, ifs),
      zn_primroot_limit : 200, zn_primroot_verbose : true;
`zn_primroot' stopped at zn_primroot_limit = 200
(%o8)                             false

Categories:  Number theory

オプション変数: zn_primroot_limit

デフォルト値: 1000

もし zn_primrootが原始根をみつけられないなら、 上限でやめます。 もしオプション変数 zn_primroot_verbose (デフォルト: false)が trueに設定されているなら、 zn_primroot_limitに到達した時メッセージが表示されます。

Categories:  Number theory

関数: zn_primroot_p  
    zn_primroot_p (x, n)  
    zn_primroot_p (x, n, [[p1, e1], …, [pk, ek]])

xが乗法群 (Z/nZ)*の原始根かチェックします。

採用したアルゴリズムは totient(n)の素因数分解を必要とします。 この素因数分解はその上時間を消費するかもしれず、いくつかの場合、 最初に素因数分解してそれから因数のリストを 3番目の引数として zn_primroot_pに渡すのが実用的かもしれません。 リストは、デフォルトオプション factors_only : falseを使って ifactors(totient(n))が返すリストと同じ形式でなければいけません。

zn_primroot, zn_order, ifactors, totientも参照して下さい。

例:

述語論理関数としての zn_primroot_p

 
(%i1) n : 14$
(%i2) units_14 : sublist(makelist(i,i,1,13), lambda([i], gcd(i, n) = 1));
(%o2)                     [1, 3, 5, 9, 11, 13]
(%i3) zn_primroot_p(13, n);
(%o3)                            false
(%i4) sublist(units_14, lambda([x], zn_primroot_p(x, n)));
(%o4)                            [3, 5]
(%i5) map(lambda([x], zn_order(x, n)), units_14);
(%o5)                      [1, 6, 6, 3, 3, 2]

オプションの三番目の引数は ifactors(totient(n))が返すリストと同じ形式でなければいけません。

 
(%i1) (p : 2^142 + 217, primep(p));
(%o1)                             true
(%i2) ifs : ifactors( totient(p) )$
(%i3) sublist(makelist(i,i,1,50), lambda([x], zn_primroot_p(x, p, ifs)));
(%o3)      [3, 12, 13, 15, 21, 24, 26, 27, 29, 33, 38, 42, 48]
(%i4) [time(%o2), time(%o3)];
(%o4)                   [[7.748484], [0.036002]]

オプション変数: zn_primroot_pretest

デフォルト値: false

もし n2, 4, p^k2*p^kなら 乗法群 (Z/nZ)*は巡回的です。 ここで pは素数で 2より大きく、 kは自然数です。

zn_primroot_pretestzn_primrootが 最小原始根を計算する前にこれらの場合の1つが起こるかどうかチェックするかどうか制御します。 zn_primroot_pretesttrueに設定されているときだけ、 これの事前テストが実行されます。

Categories:  Number theory

オプション変数: zn_primroot_verbose

デフォルト値: false

zn_primroot_limitに達したとき、 zn_primrootがメッセージを表示するかどうか制御します。

Categories:  Number theory


[ << ] [ >> ]           [Top] [Contents] [Index] [ ? ]

This document was generated by 市川雄二 on October, 5 2017 using texi2html 1.76.