Cipolla's algorithm (素数 mod 上の平方根)

はじめに

素数 mod 上の平方根を求める方法を勉強したので, 軽くまとめます

体とか同型とかそういう言葉を使わないように書いていきたいと思っています

ちゃんと勉強したい人は 37zigenさんのブログ を読んでください

やりたいこと

素数  p と 整数  a が与えられるので

 x ^ 2 \equiv a \mod p

となるような  x を一つ求めてください

存在しないならそのことを報告してください

やりかた

めんどくさいのでここから下では  \mod p をつけません (ごめんなさい)

何も言わずに  \equiv って書いてあったら  \mod p だと思ってください

とりあえず

まずは自明なものをはじきましょう

というのも, 今回使う道具は奇素数にしか対応してないです

 p = 2 のとき  x, a 0 1 のどっちかですが,

 0 ^ 2 = 0, 1 ^ 2 = 1 なので,  x = a です 簡単ですね


ついでに  a = 0 のときも先にはじきます  x = 0 ですね

ということでここから先は  p は奇素数,  1 \le a \lt p とします

平方剰余

そしたらまずは, そもそも平方根が存在するのかを判定しましょう

まず用語の説明ですが,

 x ^ 2 \equiv a となる整数  x が存在する」とき, 「  a は 法  p平方剰余 である」といいます

オイラーの規準

さて, この判定をすることができる道具として, オイラーの規準 というものがあります

これは,

 a ^ {(p-1)/2} \equiv 1 のとき  a は 法  p で 平方剰余である

 a ^ {(p-1)/2} \equiv -1 のとき  a は 法  p で 平方剰余ではない

というものです これは逆も成り立ちます


 \dfrac{p - 1}{2} がついてますね これのために p = 2 をはじいたわけです

さて,  a ^ {(p-1)/2}  1 でも  -1 でもない値になることはないのか? と思うかもしれません (僕は思いました)

これは, ならないことがわかります


フェルマーの小定理より,  p素数 a p と互いに素であるとき

 a ^ {p-1} \equiv 1 となります

ということは,

 \left( a ^ {(p-1)/2} \right) ^ 2 \equiv  a ^ {p-1} \equiv 1

ですね

二乗して  1 になる整数は  1 -1 しかないので,  a ^ {(p-1)/2} はそのどちらかになるということです

証明

 1 -1 になることはわかりましたが, なぜそれで判定できるのでしょう


まず, 「  a は 法  p で 平方剰余である 」  \Rightarrow a ^ {(p-1)/2} \equiv 1 を示します

 a は 法  p で 平方剰余なので, とある  x が存在して  x ^ 2  \equiv a です

ここで,  x p と互いに素なはずです (じゃないと  0 になってしまうので)

ということで, フェルマーの小定理より  x ^ {p-1} \equiv 1 です

 p-1 は偶数なので  x ^ {p-1} \equiv (x ^ 2) ^ {(p-1)/2} \equiv a ^ {(p-1)/2} と書けますね

ということで  a ^ {(p-1)/2} \equiv 1 です


次に,   a ^ {(p-1)/2} \equiv 1 \Rightarrow a は 法  p で 平方剰余である 」 を示します

ここで,  r p に対する原始根とします

 r は原始根なので ある自然数  k が存在して  a \equiv r ^ k と書けるはずです

よって  a ^ {(p-1)/2} \equiv r ^ {k(p-1)/2} となります

ここで  a ^ {(p-1)/2} \equiv 1 とすると,  r ^ {k(p-1)/2} \equiv 1 です

 r は原始根なので,  k(p-1)/2 p-1 の倍数である必要があります


要らない気がしますがこれも証明してみると

 k(p-1)/2 p-1 の倍数でない場合,

 k(p-1)/2 = (p-1)A + B \, ( A , B は非負整数,  0 \lt B \lt (p-1)) と書けます

ということは

  r ^ {k(p-1)/2} = r ^ {(p-1)A + B} = r ^ {(p-1)A} \cdot r ^ B

となるので

  r ^ {(p-1)A} \cdot r ^ B \equiv 1 です

ここで,  r ^ {p-1} \equiv 1 なので   r ^ B \equiv 1 であるはずですがこれは  r が原始根であることに矛盾します


ということで  k(p-1)/2 p-1 の倍数 です

つまり  k が偶数であるわけですが,  k = 2n としてみると

 a \equiv r ^ k \equiv r ^ {2n} \equiv (r ^ n) ^ 2

と書けるので  a平方根  r ^ n が存在します

証明終わり!

Cipolla's algorithm

ということで解が存在するかどうかの判定はできました

ここからは存在するものとします


ここで, いい感じの整数  b を使って

 f(x) = x ^ 2 - 2bx + a

という方程式を作ってみます

 f(x) = 0 の解は, 普通に解の公式を使うと

 x = b \pm \sqrt{b ^ 2 - a}

と書けます


さて,  b ^ 2 - a が平方剰余である場合,

ある非負整数  t が存在して  t ^ 2 \equiv b ^ 2 - a と書けるので  x \equiv b \pm t です


だからなんですか?


何もありません


何もないと困るので,  b ^ 2 - a は平方剰余じゃないとします

それはそれで困りました 平方剰余じゃないってことは  \sqrt{b ^ 2 - a} は整数じゃないです

仕方ないのでそのまま残して計算することにします


どういうことかというと,  a + b\sqrt{X}, c + d\sqrt{X} とあったとすると

 (a + b\sqrt{X}) +  (c + d\sqrt{X}) = (a + c) + (b + d)\sqrt{X}

 (a + b\sqrt{X})(c + d\sqrt{X}) = (ac + bdX) + (ad + bc)\sqrt{X}

という感じです

また, 項ごとに mod を取れることにします

 7 + 6\sqrt{X} \equiv 3 + 2\sqrt{X} \mod 4


ということで計算していきます

書くのがめんどくさいので

 b + \sqrt{b ^ 2 - a} , b - \sqrt{b ^ 2 - a} = \alpha, \beta とおきます

さて, ここでこいつらの  p 乗を考えてみます

さっきまでなら フェルマーの小定理から  a ^ p \equiv a みたいなことが言えましたが, もうそんなものは成り立ちません

 \alpha ^ p を二項定理で頑張ってみましょう

めんどくさいので  \sqrt{b ^ 2 - a} = X とおいて

 \alpha ^ p = (b + X) ^ p = \displaystyle\sum _ {k = 0} ^ p \binom{p}{k} b ^ k X ^ {p-k} と書けます


 \binom{p}{k} について考えてみます

定義より,  \binom{p}{k} = \dfrac{p!}{k! (p-k)!} と書けます

 1 \le k \le p-1 としてみると, 分母の  k! (p-k)!  p が登場しません

よって分子の  p! p が約分されないため,  \binom{p}{k} p の倍数になります

ということで mod を取ると  k = 0, p だけが残って

 \alpha ^ p \equiv b ^ p + X ^ p となります (あれ, 一年生の夢が現実になりました)


ここで,  b は普通に整数だったのでフェルマーの小定理から  b ^ p \equiv b です

 X の方は,  p が奇数なので  X ^ p = (X ^ 2) ^ {(p-1)/2} \cdot X と書けます

 X \sqrt{b ^ 2 - a} に戻すと

 (\sqrt{b ^ 2 - a}) ^ p = (b ^ 2 - a) ^ {(p-1)/2} \sqrt{b ^ 2 - a} となります


 b ^ 2 - a は平方剰余ではありませんでしたよね?

ここで, オイラーの規準を思い出してみましょう

 b ^ 2 - a は平方剰余ではないので  (b ^ 2 - a) ^  {(p-1)/2} \equiv -1 と書けます

これを代入すると

 (\sqrt{b ^ 2 - a}) ^ p  \equiv -\sqrt{b ^ 2 - a}

と書けました すっきりしましたね


そんなわけで,   \alpha ^ p \equiv b ^ p + X ^ p \equiv b - \sqrt{b ^ 2 - a}  \equiv \beta

となりました


ここで, 解と係数の関係を考えてみると,  \alpha  \beta = a です

  \alpha ^ p \equiv \beta を代入すると   \alpha ^ {p + 1} \equiv a となります


ということで,

 (\alpha ^ {(p + 1)/2}) ^ 2  \equiv a です

平方根が見つかりました!!!!


いい感じに終わった風ですが, 問題があります

ちょっと前のところを思い出してみましょう

ここで, いい感じの整数  b を使って

 f(x) = x ^ 2 - 2bx + a

という方程式を作ってみます

 b ってどうやって見つければいいんだ...


さて, この  b に求められる条件というのは

 b ^ 2 - a が平方剰余にならない

ということだけです

実は,  a が平方剰余であるとして,  bをランダムに取ってきたときに  b ^ 2 - a が平方剰余にならない確率は  \dfrac{p-1}{2p} となります (これの証明はよくわかってません...)


追記: わかりました

 b ^ 2 - a が平方剰余になるような  b がいくつあるかを数えれば, その個数を  k 個 とすると求めたい確率は  1 - \dfrac{k}{p} となります ( b の候補は  p 個あるので)

さて, 今  a は平方剰余であり  0 でないので, とある整数  m ( \not\equiv 0) を用いて  a \equiv m ^ 2 と書けます

つまり  b ^ 2 - m ^ 2 が平方剰余であるような  b の個数を数えればいいのですが,  m \not\equiv 0 なので逆元が存在して,  m ^ 2 ( (m ^ {-1} b) ^ 2 - 1 ) と書けます

 m ^ {-1} b = c とおくと,  c ^ 2 - 1 が平方剰余かどうかを判定すればいいことになります

(オイラーの規準を思い出してあげるとわかりやすいのですが, 平方剰余である数  \times 平方剰余である数 は 平方剰余である数, 平方剰余である数  \times 平方剰余でない数 は 平方剰余でない数 となります)

ここで唐突に  g(x) \equiv x ^ 2 - 2cx + 1 という関数を考えてみます  g(x) \equiv 0 の解が存在するためには  c ^ 2 - 1 が平方剰余である必要があります (上の方でやったのとだいたい同じです) これは逆も成立するので同値です

整数  t が解の一つ, つまり  g(t) \equiv 0 であるとすると,  1 が残ってしまうので  t  \not\equiv 0 であるといえます よってこいつにも逆元が存在して  c  \equiv 2 ^ {-1} (t + t ^ {-1}) と書けます

ということで  t + t ^ {-1} の値が何通りあるのか と変形できました

ここで   t + t ^ {-1} \equiv  s + s ^ {-1} とすると,   t + t ^ {-1} - s - s ^ {-1} \equiv  0 なのでこれを因数分解して  (t - s)(1 - (ts) ^ {-1}) \equiv 0 , つまり (t - s)(1 - ts) \equiv 0 と書けます

つまり  t \equiv s または  t ^ {-1} \equiv s です

したがって,   t + t ^ {-1} \equiv  s + s ^ {-1} なら  \{ t, t ^ {-1}\} = \{ s, s ^ {-1}\} となり,  \{ t, t ^ {-1}\} の個数が  k と一致します

 t p - 1 通り存在しますが, それらは  \pm 1 を除いて  t \not\equiv s \{ t, t ^ {-1}\} = \{ s, s ^ {-1}\} となるペアが存在します

ということで  k = 2 + \dfrac{p - 3}{2} = \dfrac{p + 1}{2} となります

追記終わり

ということで  \dfrac{p-1}{2p} となりますが,

これはだいたい  \dfrac{1}{2} となります

ということで  n 回連続で失敗するような確率は だいたい  \dfrac{1}{2 ^ n} くらいと言えます

何が言いたいかというと, すぐ見つかるということです

まとめ

  1.  p = 2 a = 0 の場合は early return
  2.  a ^ {(p-1)/2} を計算して解が存在するか判定
  3.  (b ^ 2 - a) ^ {(p-1)/2} \equiv -1 になるような  b を頑張って探す
  4.  (b + \sqrt{b ^ 2 - a}) ^ {(p + 1)/2} が答え

計算量は, どのステップも繰り返し二乗法を使って  \Theta(\log p) で計算できるので,

 \Theta(\log p)平方根を求めることができました

参考文献

Cipolla のアルゴリズム – 37zigenのHP

平方剰余と基本的な問題 | 高校数学の美しい物語

ルジャンドル記号とオイラーの規準 | 高校数学の美しい物語

Cipolla's algorithm - Wikipedia