元の位数を考えると群の位数計算が高速化されるが、それには高速な素因数分解が必要。「擬位数」はどの教科書にも載ってないような概念だが、ハンガリー人数学者 Babai László によって、一般の群において定義された。
p > 3 を素数、E を Fp 上の楕円曲線として、その曲線上の点の数 #E(Fp) を考える。簡単のため、以下では点の数(群の位数)を単に #E で表す。
p が小さい場合、直接的に #E を数えられる。ルジャンドル記号(実装上はヤコビ記号)の計算になるが、平方剰余の一覧表を生成して、それをヤコビ記号のキャッシュのように使ってもいいだろう。
p が大きい場合、Schoof のアルゴリズムとその発展形が使われる。
一方、p が 229 より大きく約20桁以内の場合、Shanks–Mestre のアルゴリズムを使うと、初等的な整数計算だけで手軽に #E を決定できる。その方法は baby-steps giant-steps(baby-step giant-step と呼ばれることも多い)の一種で、曲線上の(任意に選択された)点 P に対して、ハッセ区間内の整数を次々と掛けることに相当する。例外的なケースでも完全対応できるようにするためには、与えられた曲線とその2次ひねりの2種類の曲線を考える必要がある。実際には、与えられた曲線だけを考えれば十分であることが多い。
以下では、上記の文脈における baby-steps giant-steps を指して、単に BSGS と言うことにする。
BSGS にはいくつかの高速化が考えられるが、この記事では次の事柄を問題にしたい。
上記のショートカットを使った場合、仮に問題となる素因数分解を比較的素朴な方法で行ったとしても、通常、2~3倍程度の高速化が得られる。第一印象としては、この手を使った方が良いように思える。しかしながら、m の値によっては、その素因数分解を素朴な方法(試行除算)で行うと、時間がかかり過ぎる。よって、このショートカットを使うためには、m をいつでも手早く素因数分解できる仕組みを用意するか、または、「m の分解が難しいと分かったら、ショートカットをやめる」というフォールバックを準備しておく必要がある。
現実的には、次の二つを用意するのが最も良いだろう。第1は、素数性を高速に判定するアルゴリズム。これによって、素因数分解が完了したかどうかを素早く確認できるようにする。第2は、何らかの「速い」因数分解アルゴリズム。BSGS を行うときには、当然、楕円曲線の群演算を行う関数は実装されている(そうでなければ、BSGS 自体ができない)。それを再利用して、ミニECM を実装してしまえればいいだろう。…これら二つがあれば、20桁程度の数の分解は問題なくできるので、実際上、いつでも高速なショートカットが成り立つ。
要するに、ショートカットに伴う新たな問題は、微妙な落とし穴ではあるけれど、ちゃんと対応すれば比較的簡単に克服できる。
これで問題は解決なのだが、この件と関連して、ちょっとした面白い事実がある: 「多くの場合、そもそも m を完全に素因数分解する必要はない」ということだ。具体的には、位数の代わりに擬位数を使うことができる。
擬位数を使う計算は、直ちに高速化に結び付くわけではない。実際、m の素因数分解を省略してその部分の計算を節約できても、全体としては速くならない。BSGS の部分がボトルネックとなり、m の分解の手間は(常に試行除算だけを使うといった暴挙に出ない限り)もともと無視できるからだ。けれど、擬位数を使う計算は、「素因数分解が完了していないとき、完了しているふりをして、点の位数を計算すると何が起きるか」というような議論であり、それ自体として非常に面白い。
§2 では、m を素因数分解して普通に点の位数を求めて、そこから #E を導く例を提示する。§3 では、同じ例題について、m の素因数分解を完了しなくても、擬位数を用いて同様の計算が成り立つことを示す。これらの具体例の後で、擬位数に基づくアルゴリズムの一般形を検討する(§4)。
出発点として、このセクションでは「普通」の計算例を提示する。
20桁の素数 p = 86656268566282183151 を考える。これは M149 = 2149 − 1 の第1因数だが、ここでは単に任意の素数の例。
Fp 上の曲線 E: y2 = x3 + 10x + 1 の位数 #E は、ハッセの定理によって示される一定の区間
に属する。これをハッセ区間と呼ぶ。区間の幅は:
整数区間としての最小値・最大値は:
曲線 E の(勝手に選んだ)点 P = (0, 1) について、[m]P = を満たすスカラー m ∈ η を BSGS により検索すると、
が見つかる。この m はハッセ区間の中央付近にあるため、真ん中から両端に向かって探すと、比較的すぐ検出されるだろう。真ん中から探し始めるのは、実装上、そこから見たプラス側の点とマイナス側の点をペアで計算するのが効率的だから。同時に「佐藤・Tate予想」に基づく合理的な検索にもなっている(位数がある確率が高いと分かっている場所から探す)。
BSGS そのものについては、多くの文献に記されているので、ここでは説明を省略する([2], pp. 223–4 参照)。
[2] にあるような基本型の BSGS では、m がすぐ見つかっても、やっと見つかっても、どっちにしても全数検索を完了する必要があり、計算時間を短縮できない。この制約から逃れるため、m が見つかったら、P の位数計算に切り替えよう。(1) の右辺の素因数分解を利用して、P の位数 | P | を決定すると:
(1) から (2) を得るには、単に [m]P = の m の各素因数の指数を少しずつ減らしてみて、[m′]P = となる最小の m′ を決定すればいい。とりあえず計算するだけなら簡単にできる。
ここでわれわれは、暗黙のうちに (1) の素因数分解が素早くできるということを仮定している。上記の m の分解はそれほど難しくないものの、一般の m については、同じ20桁でも分解しやすい場合と分解しにくい場合がある。
さて、群の位数 #E は点の位数 | P | の整数倍で、区間 η の範囲内にある。#E の有力候補は m 自身だが、この場合 | P | > 4√p が成り立っているので、実際に #E = m = 7| P | となる。なぜなら、もし #E が m より大きければ #E ≥ 8| P | だが、これでは η の範囲外になってしまう。同様に、#E ≤ 6| P | もあり得ない。別の言い方をすると、#E は | P | の倍数だが、この | P | は 4√p より大きいので、| P | の倍数のうちハッセ区間内に入るものは、たかだか1個しかない(群論的な性質とハッセの定理により、必ず1個は存在する)。
上記の「元の位数の計算法」は、米国の Sutherland (サザランド)が高速位数アルゴリズムと呼んでいるもので、一般の群において成立する。#E(Fp) の計算への応用は、2015年の [7], Algorithm 8.2 において講義の題材となっており、[8], Problem 3 は関連する演習問題である。一般の群の元に対する高速位数アルゴリズムは、2007年の [5], Chapter 7 において、検討されている。
BSGS と高速位数アルゴリズムを組み合わせるというアイデアは Sutherland の発見ではなく、1990年代に Henri Cohen (アンリ・コアン)が [1], Algorithm 7.4.12 として記述している。Cohen()はフランスの数学者で、PARI/GP()の開発者。注意点として、1993年の初版では定理の成立条件が不正確であり、p > 19 となっている。正誤表によると、第2版以降では「定理としては p > 457 だが、アルゴリズムを少し変えれば p をかなり小さくできる」と改訂された。
実際、「曲線とその2次ひねりの少なくとも一方に | P | > 4√p となる点がある」という保証を求めた場合、最善の評価は p > 457 である(Mestre の定理)。アルゴリズムの目的上は、| P | > 4√p の代わりに m ± | P | ∉ η とすれば(言い換えれば m − | P | < ηmin かつ m + | P | > ηmax とすれば)十分で、その場合、条件が p > 229 に緩和される(Mestre–Schoof の定理)。
Cohen が初版で書いた条件は、必ずしも完全な誤りではない。[6] によると、定理の内容を少し変えれば p > 29 にでき、少数の例外ケースを場合分けすれば、もっと小さくもできるようだ。
§2 と同じことを「普通」でない方法でやってみる。
上記の (1) の計算において、素朴な試行除算で素因数分解を行うとすると、約1000万までの全部の素数で割ってみることが必要になる。BSGS 自体の手間に比べればささいな労力かもしれないが、数が大きくなると試行除算のコストは急激に増大する。
直球勝負としては、「p − 1」法などの強力な分解手段を使えば良さそうだ。実装の手間という点では、とりあえず試行除算してみて、それで簡単に分解できない場合には、ECM を使うのが手軽そうだ。しかしここでは逆に、素因数分解の「手抜き」を考える。話を具体的にするため、220 = 1048576 以下の素数による試行除算を行い、その範囲で分解完了できない場合、余因数を放置する。
そうすると (1) の分解は完了せず、次のように余因数 c = 80108370572519 = 8506843 × 9416933 が残る。
20桁程度の数を試行除算で素因数分解する場合、分解が完了したかどうかは、余因数に対する Miller–Rabin テストで簡単に判定できる。2015年の発見 [9] によると、最初の12個の素数を底として最大12重のテストを行えば、23桁以内の数については確定的な(擬素数フリーの)素数判定ができる。同様に、最初の13個の素数を使えば、25桁の 3.317×1024 までは確定的に判定できる。詳細については OEIS: A014233 を参照。25桁より大きい場合の判定については、今回は考えない。
c は合成数だが、それを含めて、上記の各因数は互いに素である。ハンガリー出身の数学者 Babai László (ババイ・ラースロー(): Babai が姓)は1999年に出版された [4] の中で、このような
をみなし素因数(pretend-primes)と呼んでいる。「互いに素な因数」と言い換えてもいいだろう。
これを基に最初と同様の計算を行うと、次が得られる:
すなわち、因数 7 をなくしても依然
となるが、それ以外の因数は1乗たりとも減らせず、減らせば積は にならない。例えば:
このことから、| P | の正確な値は、2 × 77267 × c′ の形を持つことが分かる。ここで c′ ≠ 1 は c の約数であり、c の素因数1個以上から成る。ところが、仮定により c には 220 = 1048576 以下の素因数がない。ゆえに:
この r1 は 4√p より大きいので、最初の例と同様に、#E = m が示される。m の素因数分解は完了しなかったし、従って点の正確な位数も確定できなかったけれど、点の位数の下界が分かり、目的を達することができた。
#E = m ではない可能性があるためには、
であることが必要。この条件が成り立てば、m ± | P | の少なくとも一方 m2 も #E の候補となり得る。実際、m2 は | P | の倍数だから、ハッセ区間内にある限りにおいて、それが群の位数である可能性がある。逆に、上記の不等式が成り立たないなら、ハッセ区間に | P | の倍数は1個しかない。
ハッセ区間内に入るか・入らないかという部分を具体的に考えると、より精密な評価が得られる:
これが成り立ってしまうと、#E = m と断定できない。逆に
であれば、#E = m が確定する。c′ が c の因数であることに注意すると、最後の不等式を成り立たせるためには、「c には、右辺の値以下の小さい素因数はない」ということを示せばいい。一般的に言えば、2 × 77267 は m′/c であり、「c には、次の B2 以下の素因数はない」ということを示せば、#E = m となる。
ここで考えている例においては、
この例において、c には B1 = 1048576 以下の素因数はないのだから、B2 以下の素因数がないことは明らか。いつもこんなに簡単に解決するとは限らないが、この例では、うまくいった。
Babai は、上記の m′ のようなものを P の擬位数(pseudo-order)と呼んでいる。擬位数は、
であり、その点の(より一般的に言えば、その元の)位数の倍数である。Sutherland も、この擬位数という用語を使っている([5], p. 109)。
楕円曲線の位数に関する上記の計算例は、場当たり的な奇策のように見えるかもしれない。しかし Babai たちは、一般の群について、次のように指摘している([4], p. 15)。
多項式時間では、元の位数を決定することすら、できないかもしれない。位数の決定は、整数の因数分解という計算数論の別の難問(多項式時間では解けないと信じられている)と結び付いているからだ。(中略)しかしながら、大抵の目的においては、元の正確な位数を使うことを回避し、みなし素因数の集合に基づく元の擬位数で間に合わせることができる。
Babai の発想の萌芽は、1997年の [3] に見られる。
筆者は Babai の論文を読んで上記の計算を行ったわけではなく、最初は実際、場当たり的にこの方法を思い付いた。後から考えてみると、擬位数を経由する #E の計算は、Babai が指摘したことの実践例といえるだろう。
Fp 上の楕円曲線 E の点 P について、ハッセ区間を探索して [m]P = となる m が見つかったとする。
擬位数を使う場合、m の素因数分解は不完全でも構わない。とりあえず、ある定数 B1 以下の素因数だけを考えればいい。もし m の素因数のうち、最大素因数以外の全てが B1 以下なら、上記の範囲を考えるだけで完全な素因数分解が成り立ち、以降は普通のアルゴリズムとなる(最大素因数は、余因数として自動的に検出される)。以下では、そうならない場合、つまり m に含まれる B1 以下の素因数を全て検出したとき、余因数が依然として合成数である場合を考えよう。例えば、次のように進む:
m に含まれる B1 以下の素因数を pi (i = 1, 2, …, k) とし、余因数を c とせよ。仮定により、c は B1-ラフである(B1 より大きい素因数だけを含む)。次の形の「不完全素因数分解」が成り立つ(ui ≥ 1 は整数):
各 pi と c は互いに素であり、これらは「みなし素因数」を成す。
P の擬位数 m′ は、次の形を持つ(0 ≤ vi ≤ ui, 0 ≤ w ≤ 1 は整数):
これに対応するもし w = 0 なら m′ = | P | であり、通常のアルゴリズムが使えるので、w = 1 と仮定する。擬位数の不完全素因数分解における余因数は、引き続き c1 = c である。
B2 = ⌊c max(m − ηmin, ηmax − m)/m′⌋ と置け。もし B2 ≤ B1 なら m = #E で終了。素因数分解が未完了のままで、計算を完了させたことになる。
B1 < B2 なら、素因数分解を再開しなければならない(ただし困難があれば、いつでも諦めて普通の BSGS に戻ることもできる。現実的には、ここで ECM に切り替えるべきだろう)。(B1, B2] にある各素数で、小さい順に c を割ってみる。もしどれでも割り切れないなら、m = #E で終了。これは面白いパターンであり、「素因数分解を再開したが、新たな素因数は見つからなかった。にもかかわらず、問題自体は解決した」という形になっている。
c の素因数 q が見つかった場合には、そこでいったん停止して、余因数 c/qu を検討する(ここで u は整数で、c が q で何回割り切れるかを表す。ほとんどの場合 u = 1)。もし c/qu が合成数でないなら、m は完全に素因数分解されたことになり、通常のアルゴリズムが使える。この場合、この手順に従うと「最初から真面目に完全分解した方が早かった。手抜きをしたら二度手間になってしまった」ということになる。
一方、途中でc/qu が依然として合成数になっている可能性が考えられる。その場合、素因数 pk+1uk+1 = qu が増えた状態で (1) に戻ることになる。B1 := q, c := c/qu, k := k + 1 となる。
もう一つのパターンとして、少なくとも理論上は、余因数引き続き p = 86656268566282183151 として、Fp 上の曲線 E: y2 = x3 + ax + 1 と、その曲線上の点 P = (0, 1) を考える。ハッセ区間の下限と上限は:
BSGS により、ハッセ区間内を検索して [m]P = を満たす m を発見し、それを分解する。このとき B1 = 220 = 1048576 以下の素因数だけを考える。m が B1 以下の素因数と素数の余因数しか持たないなら、完全な素因数分解になるので、そうでない場合について考察する。a ∈ [1, 100] の範囲には、そのような例が12種類ある。それらについて、高速位数アルゴリズムによって擬位数 m′ を求め、前述の公式によって B2 を定義する。
a = 10 の場合は §3 で取り上げたが、
最後の行の計算は、⌊r/(2 × 77267)⌋ としても、同じことである。
この B2 は B1 以下なので、m が群の位数であることが(言い換えれば、この楕円曲線の上には m 個の点があることが)直ちに確定される。
同様に、a = 17 の場合:
やはり B2 は B1 以下であり、同じ結論になる。
a = 25, 26, 76 の場合も、全く同様。ここでは、擬位数が m に等しくなる例として、a = 26 を取り上げる:
a = 81 の場合:
B1 < B2 となったので、(B1, B2] = (1048576, 1274486] 間の各素数で c を割る。どれでも割り切れず m の素因数分解は進展しないが、B1 = 1274486 と再設定したのと同じことになるため、この時点で B2 ≤ B1 が満たされ、m が群の位数であると結論される。
a = 95, 96 の場合も同様に進むが、それぞれ B2 = 7481965, 64514143 となり、後者では約6500万までの試行除算が必要とされる。c に隠れている実際の素因数は約1~2億なので、完全な素因数分解よりはいくぶん計算量が節約される。
実際には、このような場合(特に B2 が大きい場合)、m を分解するための強力なアルゴリズムに切り替えるか、または m の分解を諦めて、BSGS に戻った方がいいだろう。
a = 60 の場合:
B2 は約3600万となるが、(B1, B2] 間の素数で順に c を割ると、300万台の素因数 q = 3834557 が新たに発見され、約300億の余因数 c/q も素数となる。擬位数を使ったことで結果的には遠回りになってしまったが、完全な素因数分解が成立した。この場合、正確な | P | を求めることも容易だが、
さえ確認すれば、| P | に約300億という大きさの素因数が含まれることが分かり、それだけでも m = #E と結論できる。
a = 62 も同様に進むが、B2 = 468816154 となり、素朴にやると約5億までの素数表が必要となる。実際には、わずか100万台で因数 q = 1798183 が見つかり、素因数分解が完了する。
a = 73, 74 については、B2 はそれぞれ約2000万、1億となるが、実際にはそれぞれ約1000万、200万の素因数が見つかって、分解完了する。
これらの例では、擬位数を使う作戦が奏功せず、結局、素因数分解を完了させることになった。このような場合、一般的に言って、m を分解するための強力なアルゴリズムに切り替えるか、または m の分解を諦めて、BSGS に戻った方がいいだろう。
例えば、ECM に切り替えたとして、その場合も、新しい素因数の発見は必要でも必ずしも完全な素因数分解は必要なく、擬位数計算で十分である可能性がある。
前回の「楕円曲線で因数分解」の §4 などでは、PARI/GP を使って位数計算を行った。その後、自力で位数計算を行えるように とりあえず BSGS を実装し、Sutherland に従って「高速位数アルゴリズム」と組み合わせるショートカットを試みた。このとき「高速位数アルゴリズムのための素因数分解が、面倒なこともある」という問題に直面した。
Sutherland も Henri Cohen も、「ショートカットあり」の BSGS を紹介しつつ、ショートカット部分の因数分解の手間を明示的に評価していない。この点は、実装上・学習上の落とし穴となり得る。
メールで質問してみたところ、Henri Cohen 先生から「実用上、普通は素朴な因数分解でも十分だが、確かに素朴な方法で因数分解すると、そのコストは BSGS の複雑性 O(p1/4) を上回ってしまう」という明快な回答をいただいた。「解決法としては、ECM を実装するのが恐らく最も簡単だろう」ということで、筆者の考えていた選択肢と同じだった。Sutherland 先生からも実質的に同内容の返信をいただいた。ニュアンスの違いとして、Henri Cohen は「厳密に言うと、自分のアルゴリズムにはそういう意味では穴がある」という感じだが、Sutherland は「数学的に間違ったことは言っていないから、あのアルゴリズムに穴はない。素因数分解の問題は実装の問題であり、アルゴリズムの欠陥ではない」という感じだった。講義ノートには「この素因数分解は比較的簡単」と書いてありながら、実際には「単純な試行除算では駄目」という制約がある。この点は不透明ではあるが、「比較的簡単な方法が存在する」という主張は数学的に正しい。講義ノートは教科書ではないので、技術的細部が省略されている場合もあるだろう。
「完全分解しなくても、いけるのでは」というのは楽しい再発見だった。§5 で考えた100種類の楕円曲線のうち88種類では、最初から真の位数が得られた。それ以外の12種類のうち5種類では、擬位数を使うことで、完全な素因数分解を容易に回避できた。ただし、常に素因数分解を回避できるわけではなく、このトリックは上記の「落とし穴」の解決法とはならない。
BSGS の計算量が大きいので、素因数分解の部分で計算量を減らせたとしても、トータルの計算時間にはほとんど影響しない。そもそも、もっと強力なアルゴリズムが既にあるのだから、大きな数に対する BSGS を高速化することそれ自体、あまり意味がないのかもしれない。けれど、実用性の大小にかかわらず、擬位数計算というアプローチは、それ自体として面白いと感じられた。