Tonelli-Shanks アルゴリズムの理解と実装(1)

はじめに

奇素数$p$と$0$でない整数$n$が与えられた際に

$$ x^2 \equiv n\ {\rm mod}\ p$$

なる$x\not \equiv 0$が存在するならば、$n$は$p$を法として**平方剰余(Quadratic Residue)であると言い、そのような$x$が存在しないならば$n$は平方非剰余(Quadratic Nonresidue)**であると言います。

例えば、$n=2$は$p=7$を法として平方剰余($x=3$などが解)ですが、$n=3$は平方非剰余です。

さて、平方剰余かどうかの判定はそれほど難しくありません。

まずは、Legendre(ルジャンドル)記号を以下で定義します。

$$\left(\begin{array}
nn\
p
\end{array}\right) := n^{\frac{p-1}{2}}$$

この時、以下が成立します(証明は少し難しいのでここでは省略します)。

$$\left(\begin{array}
nn\
p
\end{array}\right) \equiv \begin{cases}
1\ {\rm mod}\ p\Leftrightarrow n は平方剰余\
-1\ {\rm mod}\ p \Leftrightarrow n は平方非剰余
\end{cases}$$

実際、$n = 2, p = 7$なら

$$\left(\begin{array}
nn\
p
\end{array}\right) = 2^{\frac{7-1}{2}} = 8 \equiv 1\ {\rm mod}\ 7$$

であり、$n = 3, p = 7$なら

$$\left(\begin{array}
nn\
p
\end{array}\right) = 3^{\frac{7-1}{2}} = 27 \equiv -1\ {\rm mod}\ 7$$

です。

この記事では、$n$が平方剰余であると分かった場合に、その平方根$x$を求める方法について書いていきます。以下、特に断りがなければ$\equiv$は全て${\rm mod}\ p$です。

pが4で割って3余る素数の場合

この場合は実は簡単で、

$$x = \pm n^{\frac{p+1}{4}}$$

が答えです。実際、$n$は平方剰余であるという前提なので

$$\left(\begin{array}
nn\
p
\end{array}\right) = n^{\frac{p-1}{2}} \equiv 1$$

が成立し、

$$x^2 = n^{\frac{p+1}{2}} = n^{\frac{p-1}{2}} \cdot n \equiv n$$

となります。例えば、$n = 2, p = 7$なら

$$x = \pm 2^{\frac{7+1}{4}} = \pm 4 \equiv 3, 4\ {\rm mod}\ 7$$

となります。

pが4で割って1余る素数の場合

ここが今回取り上げるTonelli-Shanksのアルゴリズムになります。

記号の命名は基本的にWikipediaに準拠しています。

Step 1.

$p-1$を$2$で割れるまで割り算し、

$$p-1 = Q \cdot 2^S$$

という形に変形します($Q$は奇数で、$S$は正の整数です)。

Step 2.

$p$を法として平方非剰余である数をランダムに選び、それを$z$とします。

「ランダムに」と言っていますが、$1$から$p-1$までの間に平方剰余と平方非剰余の個数は半分ずつなので¹、$z = 1, 2, …$と小さい順に総当たりしても計算量的には問題ありません。

なお、平方非剰余かどうかは上記のLegendre記号を使えば確かめられますし、当然その定義上、

$$z^{\frac{p-1}{2}} \equiv -1$$

です。

Step 3.

まず4つの数$M_0, c_0, t_0, R_0$を以下のように定義します。

$$\begin{cases}
M_0 = S\
c_0 = z^Q\
t_0 = n^Q\
R_0 = n^{\frac{Q+1}{2}}
\end{cases}$$

Step 4.

以下のようにループ文および漸化式を立てます。

1. もし$t_i \equiv 1$なら、$\pm R_i$が$n$の平方根であり、ループ文を抜けて終了。

2. そうでない場合は、以下のように値を更新する。

$$\begin{cases}
M_{i+1} = \left(\left(t_i\right)^{2^{j}}\equiv 1 を満たす最小のj, ただし0 < j < M_i\right)\
\
c_{i+1} = \left(c_i\right)^{2^{M_i - M_{i+1}}}\
\
t_{i+1} = t_i \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}}\
\
R_{i+1} = R_i \cdot \left(c_i\right)^{2^{M_i - M_{i+1}-1}}
\end{cases}$$

…一気に見にくくなりました。

• 「もし$t_i \equiv 1$なら、$\pm R_i$が$n$の平方根であり」とはどういう事か？

• 無限ループになることはないのか？

• 実際にコーディングして動かしたい / 例を見せてくれ

などなどあると思いますが、順番に解説していきます。

準備

漸化式をいくつも立てていますが、実は$i$の値によらず以下が成立します。

$$\begin{cases}
\left(c_i\right)^{2^{M_i -1}} \equiv -1\
\
\left(t_i\right)^{2^{M_i -1}} \equiv 1\
\
\left(R_i\right)^2 \equiv t_i \cdot n
\end{cases}$$

順番に数学的帰納法を使って確認していきます。

1つ目の式

$$\left(c_i\right)^{2^{M_i -1}} \equiv -1$$

を示します。まず、$i = 0$ならば、

$$\left(c_0\right)^{2^{M_0 -1}} = \left(z^Q\right)^{2^{S -1}} = z^{Q \cdot 2^{S -1}}$$

ここで、Step 1.より$\frac{p-1}{2} = Q \cdot 2^{S-1}$なので、

$$z^{Q \cdot 2^{S -1}} = z^\frac{p-1}{2}$$

これはStep 2.で定義したので$-1$に等しくなります。

つぎに、$i$で成立するとして$i+1$で成立することを示します。

$$\left(c_{i+1}\right)^{2^{M_{i+1} -1}} = \left(\left(c_i\right)^{2^{M_i - M_{i+1}}}\right)^{2^{M_{i+1} -1}} = \left(c_i\right)^{2^{M_i - 1}}$$

となり、$i$で成立すると仮定しているのでこれも$-1$になります。

2つ目の式

$$\left(t_i\right)^{2^{M_i -1}} \equiv 1$$

を示します。$i = 0$ですが、1つ目の時と同様にして

$$\left(t_0\right)^{2^{M_0 -1}} = \left(n^Q\right)^{2^{S-1}} = n^{Q \cdot 2^{S-1}} = n^\frac{p-1}{2}$$

$n$は平方剰余であるという仮定なので、これはLegendre記号の定義から$1$に等しくなります。

次に、$i+1$についてですが、これは少し込み入った証明となります。

$$\left(t_{i+1}\right)^{2^{M_{i+1} -1}} = \left(t_i \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}}\right)^{2^{M_{i+1} -1}} = \left(t_{i}\right)^{2^{M_{i+1} -1}} \cdot \left(c_{i}\right)^{2^{M_i -1}}$$

ここで、$\left(t_{i}\right)^{2^{M_{i+1} -1}}$に注目します。

これを二乗すると、$\left(t_{i}\right)^{2^{M_{i+1}}}$になり、Step 4.の$M_{i+1}$の定義から$1$に等しくなります。

では、二乗すると$1$になる数は？と言われれば、もちろん$1, -1$の$2$つだけです²。

しかし、今回$\left(t_{i}\right)^{2^{M_{i+1} -1}} \equiv 1$…($\ast$)とはなりません。なぜなら、もしこれが成立してしまうと、「Step 4.で

$$\left(t_i\right)^{2^{j}}\equiv 1 を満たす最小のj$$

を$M_{i+1}$と定義したにもかかわらず、$j = M_{i+1}-1$でも上の式($\ast$)が成立してしまっている」という矛盾が生じてしまうためです。

よって、$\left(t_{i}\right)^{2^{M_{i+1} -1}} \equiv -1$となります。

また、$\left(c_{i}\right)^{2^{M_i -1}}$は$1$つ目の式から$-1$と証明したばかりなので、結果として、

$$\left(t_{i}\right)^{2^{M_{i+1} -1}} \cdot \left(c_{i}\right)^{2^{M_i -1}} \equiv (-1) \cdot (-1) = 1$$

となります。

3つ目の式

$$\left(R_i\right)^2 \equiv t_i \cdot n$$

を示します。$i = 0$ならば、

$$\left(R_0\right)^2 = \left(n^{\frac{Q+1}{2}}\right)^2 = n^{Q+1} = n^Q\cdot n = t_0 \cdot n$$

また、$i+1$の場合は

$$\left(R_{i+1}\right)^2 = \left(R_i \cdot \left(c_i\right)^{2^{M_i - M_{i+1}-1}}\right)^2 = \left(R_i\right)^2 \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}}$$

$i$で成立、つまり$\left(R_i\right)^2 \equiv t_i \cdot n$を仮定しているので

$$\left(R_i\right)^2 \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}} \equiv t_i \cdot n \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}}$$

定義より$t_{i+1} = t_i \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}}$であるので

$$t_i \cdot n \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}} = t_{i+1} \cdot n$$

以上で$3$つの式全てを証明できました。

各疑問点の解消

まず

「もし$t_i \equiv 1$なら、$\pm R_i$が$n$の平方根であり」とはどういう事か？

から説明します。とは言っても、上で示してきた式の$3$つ目

$$\left(R_i\right)^2 \equiv t_i \cdot n$$

を見れば明らかです。つまり、$t_i$を更新していき、どこかのタイミングで$t_i\equiv 1$になれば、その時点で上の式は

$$\left(R_i\right)^2 \equiv n$$

と同じ意味になります。

次に、

無限ループになることはないのか？

です。つまり、$t_i = 1$になるような$i$はないのではないか、ということです。

それにはまず、$M_i$の単調減少性について説明する必要があります。$M_{i+1}$を

$\left(t_i\right)^{2^{j}}\equiv 1 を満たす最小のj, ただし0 < j < M_i$

と定義していましたが、そもそも本当に$0 < j < M_i$にそのような$j$が存在するのでしょうか？

実は、その点については、上で示してきた式のうち$2$つ目

$$\left(t_i\right)^{2^{M_i -1}} \equiv 1$$

が保証してくれています。つまり、$j = 1, 2, …$と調べていき、運が悪くてずっと$\left(t_i\right)^{2^{j}}\not \equiv 1$だったとしても、$j = M_i -1$の時には$\left(t_i\right)^{2^{j}}\equiv 1$になり$M_{i+1}$としての条件を満たします。

よって、$M_{i+1}$は必ず$M_i -1$以下となり(添え字があるのでわかりにくいですね…)、

$$M_{i+1} \leq M_i -1 < M_i$$

これで、$M_i$の単調減少性が言え、いつかは$M_{end} = 1$となる時が来ます。

この時、上で示してきた式のうち$2$つ目

$$\left(t_i\right)^{2^{M_i -1}} \equiv 1$$

を用いると、

$$\left(t_{end}\right)^{2^{M_{end} -1}} \equiv 1 \Rightarrow t_{end} \equiv 1$$

$t$の値が$1$になったので、これで無限ループすることはなくなりました。

例

Wikipediaの例をそのまま引用します。

$x^2 \equiv 5 \ {\rm mod}\ 41$

なる$x$を求めます($n = 5, p = 41$)。まず、Legendre記号を使って

$$\left(\begin{array}
55\
41
\end{array}\right) = 5^{\frac{41-1}{2}} = 5^{20} \equiv 1\ {\rm mod}\ 41$$

と、$5$が平方剰余であることを念のため確認しておきましょう。

Step 1.

$p-1 = Q \cdot 2^S$、つまり$41-1 = 5 \cdot 2^3$とします。$Q = 5, S = 3$です。

Step 2.

平方非剰余である数$z$を選ぶ。Legendre記号を使うと、$1, 2$は平方剰余であるが、$3$は

$$\left(\begin{array}
33\
41
\end{array}\right) = 3^{\frac{41-1}{2}} = 3^{20} \equiv -1\ {\rm mod}\ 41$$

と、平方非剰余であるので$z=3$とする。

Step 3.

$$\begin{cases}
M_0 = S\
c_0 = z^Q\
t_0 = n^Q\
R_0 = n^{\frac{Q+1}{2}}
\end{cases}$$

つまり、

$$\begin{cases}
M_0 = 3\
c_0 = 3^5 \equiv 38\
t_0 = 5^5 \equiv 9\
R_0 = 5^{\frac{5+1}{2}} \equiv 2
\end{cases}$$

Step 4.

ここから、$t_i=1$になるまでループに入ります。

$$\begin{cases}
M_{i+1} = \left(\left(t_i\right)^{2^{j}}\equiv 1 を満たす最小のj, ただし0 < j < M_i\right)\
\
c_{i+1} = \left(c_i\right)^{2^{M_i - M_{i+1}}}\
\
t_{i+1} = t_i \cdot \left(c_i\right)^{2^{M_i - M_{i+1}}}\
\
R_{i+1} = R_i \cdot \left(c_i\right)^{2^{M_i - M_{i+1}-1}}
\end{cases}$$

$i=0$の時

$$\begin{cases}
M_{1} = \left(\left(9\right)^{2^{j}}\equiv 1 を満たす最小のj\right) = 2\
\
c_{1} = \left(38\right)^{2^{3 - 2}} \equiv 9\
\
t_{1} = 9 \cdot \left(38\right)^{2^{3 - 2}} \equiv 40\
\
R_{1} = 2 \cdot \left(38\right)^{2^{3 - 2-1}} \equiv 35
\end{cases}$$

$i=1$の時

$$\begin{cases}
M_{2} = \left(\left(40\right)^{2^{j}}\equiv 1 を満たす最小のj\right) = 1\
\
c_{2} = \left(9\right)^{2^{2 - 1}} \equiv 40\
\
t_{2} = 40 \cdot \left(9\right)^{2^{2 - 1}} \equiv 1\
\
R_{2} = 35 \cdot \left(9\right)^{2^{2 - 1-1}} \equiv 28
\end{cases}$$

$t_2 = 1$なので、答えは$\pm R_2 = 28, 13$となります。

実際、確かに$13^2 \equiv 28^2 \equiv 5 \ {\rm mod} \ 41$は成立します。

実装

長くなってしまったので、こちらに回します。

1. 簡単な証明を載せておくと、$i^2 \equiv (p-i)^2$なので、$1$から$p-1$までの数を二乗した平方数のバリエーションとしては多くても$\frac{p-1}{2}$個です。しかし、「多くても」と言っていますがこれ以外に重複はありません。これは$i^2 \equiv (i+k)^2$を仮定して背理法で示せます。 ↩

2. 素数$p$でしか成立しません。合成数なら、$x^2 \equiv 1\ {\rm mod}\ 24$なる$x$として、$x = 1$と$x = -1 \equiv 23$以外にも$x=5$などが考えられます。素数$p$で$x = \pm 1$しか解をもたない理由は、$x^2-1$を因数分解すると容易に証明できます。 ↩