二次篩法 (QS)


改めて GNFS の勉強をし直してみたくなったので脳内で分かってるつもりになっている事を確認しつつ実装していく予定で書き始めてみました。まずは理論も実装も簡単な QS からスタートです。

篩法

概要

まずは今後やっていく篩系の素因数分解アルゴリズムに共通する概要を説明します。各篩法では最終的に因数分解したい合成数 N に対して

x^2 \equiv y^2 \pmod{N}

となる非自明な整数組 (x,y) を見つけることになります。
すると x^2-y^2=(x+y)(x-y)\equiv 0 から x\pm yN の素因数をバラバラに持っていることを期待して x\pm yN の最大公約数を求めます。その結果が 1N になった場合は残念ながら分解失敗となりますが、他の値が出ればそれは N の(素)因数ということで成功となります。

例 (N=187)

例えば N=187 を例として以下の3式を見てみましょう。(以下\pmod{187}は省略)

\begin{align*} 4^2 & \equiv -171 = (-1) \cdot 3^2 \cdot 19\\ 5^2 & \equiv -162 = (-1) \cdot 2 \cdot 3^4\\ 15^2 & \equiv 38 = 2 \cdot 19 \end{align*}

すると右辺は素因数分解の結果から全部かけると平方数になることが分かるので

\begin{align*} (4\cdot5\cdot15)^2 &\equiv ((-1) \cdot 2 \cdot 3^3 \cdot 19)^2\\ 300^2 &\equiv 1026^2 \end{align*}\\ (x,y)=(300,1026)

という計算を経て {\rm GCD}(N,x+y)=17 と分解できました。

では一般の合成数 N に対してこの 3 式のようなものはどのようにすれば効率的に見つけられるのでしょうか、というところが篩法の突き詰めていく部分になります。

二次篩法

概要

二次篩法 (QS; Quadratic Sieve) は以下のように大きく4つのステップから構成されます。

  1. パラメータの選択

  2. 有効式の選択

  3. 分解

  4. の篩ステップでは 1. で決めた範囲 [a, b) にある整数 x_k について

f(x) = x^2 - N

という二次式の値の素因数分解を考えます。ここで [a,b) として \sqrt{N} を含む範囲を設定しているとすると、f(x_k)O(\sqrt{N}) の大きさになるのでそこそこ小さな素数だけで分解できる可能性があります。

上の N=187 について [a,b)=[4,16) として見てみましょう。20 以下の素数だけで分解できた数だけ書いておくと

\begin{align*} 4^2 & \equiv -171 = (-1) \cdot 3^2 \cdot 19\\ 5^2 & \equiv -162 = (-1) \cdot 2 \cdot 3^4\\ 11^2 & \equiv -66 = (-1) \cdot 2 \cdot 3 \cdot 11\\ 13^2 & \equiv -18 = (-1) \cdot 2 \cdot 3^2\\ 14^2 & \equiv 9 = 3^2\\ 15^2 & \equiv 38 = 2 \cdot 19 \end{align*}

となります。

次にステップ 3. の有効式の選択では上記のようにある程度の素数で分解できた f(x_k) の情報からかけ合わせて平方数になる組み合わせを選択します。上記の例では \bm{x}=\{4,5,15\} だけでしたが、\{5,13\}\{5,13,14\} といった組み合わせも考えられます。

そして 4. の分解ステップで既に篩法の概要で述べたように {\rm GCD}(N, x\pm y) の計算をすることで非自明な因数を見つけられることが期待されます。

アルゴリズム

ここまでは数式の上で原理を理解するための概要説明であり、具体的な手順として不明なところや素直に実行しては効率的ではないところがあったかと思いますので順に説明しましょう。

篩ステップでは f(x_k) の素因数分解を行うことになりますが、これを直接計算すると結局試し割りの連続となるので効率がよくありません。

for (Integer x = a; x < b; ++x) {
  Integer fx = f(x);
  for (int p : factor_base) {
    while (fx % p == 0) {
      exponent[x][p]++;
      fx /= p;
    }
  }
  cached_f[x] = fx;  // == 1 なら FB-smooth
}

なのでこの f(x_k) の因数分解を纏めて効率的に行います。

素数 p がある x=x_k について f(x) を割り切れたとしましょう。

p | x_k^2-N

すると x=x_k+p についても p で割り切ることができるはずです。

f(x_k+p) = x_k^2 + 2px_k + p^2 - N = f(x_k) + p(2x_k+p) \\ \therefore p | f(x_k+p)

なので f(x)p で割り切れる x を 1 つ見つけるとその前後 p 毎に p で割り切れることがわかります。また、そのような x は方程式

x^2-N \equiv 0 \pmod{p}

を満たすので、\bmod{p} における N の平方根を探ることになります。一旦詳細は省きますが、いくつかの前提のもとこのようなアルゴリズム で求めることができます。

for (int p : factor_base) {
  for (int r : ModSqrts(n % p, p)) {
    for (Integer x = a - (a % p) + r; x < b; x += p) {
      while (cached_f[x] % p == 0) {
        exponent[x][p]++;
        cached_f[x] /= p;
      }
    }
  }
}

ModSqrts(a, p)\bmod{p} での N の平方根2つを返す関数です。
さて、この方法を改めてみてみるとエラトステネスの篩と同じような作業をしていることが分かるでしょう。これが「篩系」アルゴリズムと言われる理由です。