Hatena::ブログ(Diary)

oupoの日記

2012-11-22

線形合同法であるseedが0からいくつ進めたものかを得る

12:23

sid-searchを作り直してみた - oupoの日記でも使ったアルゴリズムを解説する。

問題

線形合同法の漸化式に使われる関数(x * A + B) mod 2^Nをf(x)とおき、この線形合同法の数列は最大周期であるとする。

<a_n>を初項を0としたこの線形合同法の数列とする。すなわちa_0 = 0, a_{n+1} = f(a_n)

与えられたseed、xに対してf^i(0) = xとなるi ∈ 0...2^Nを求めたい。

思考

a_nは偶数奇数を繰り返すことからxの偶奇からiの偶奇がわかる。iの第0ビットがわかったので、同じことを繰り返して第1ビット、第2ビット...と次々に求めていけないだろうかと考える。

  • iが偶数の場合、0にf^2を何回適用すればxになるかを求めて、それを2倍すれば0にfを何回適用すればxになるかが得られる。
  • iが奇数の場合、0にf^2を何回適用すればf(x)になるかを求めて、それを2倍して1引けば0にfを何回適用すればxになるかが得られる。

この考えを使って繰り返せばよい。

f:id:oupo:20121122120936p:image

ここで、f^i(0) = xとなるiの偶奇はxの偶奇から得られたが、fの代わりにf^2を当てはめた場合のiの偶奇は得られるか?

xの第1ビットから分かるのではないかと予想できる。

実際:

na_na_n mod 4
00x000000000b00
20x6c0789660b10
40xdbffe6dc0b00
60x992415e20b10
80x895277f80b00
100x014d329e0b10
120x75ab8f540b00
140xe0d45b9a0b10
160xe69948f00b00
180x964b4cd60b10

一般に、g=f^{2^k}とおいたとき、g^i(0) = xとなるiを2で割った余りはxの第kビットに一致する。

これは<a_n>の下位kビットだけ注目したとき周期が2^kであることから証明できる。

k=1を例にとって解説してみる。

<a_n mod 4>の周期は4で、この表を考える。

na_n mod 2a_n mod 4
00b00b00
10b1
20b0
30b1

a_n mod 2の値からa_n mod 4の第0ビットは分かる。

na_n mod 2a_n mod 4
00b00b00
10b10b?1
20b00b?0
30b10b?1

a_2 mod 4が0b00だとすると周期が4なことに反してしまうのでa_2 mod 4 = 0b10である。

na_n mod 2a_n mod 4
00b00b00
10b10b?1
20b00b10
30b10b?1

よって0を初項、f^2を漸化式とする数列はmod 4で見たとき0b00, 0b10を繰り返す。

プログラム

# f^m(0) = sとなるm ∈ 0...2^Nを求める
def seed_to_index(s)
  seed_to_index0(s, 0)
end

# g = f^{2^n}としてg^m(0) = sとなるm ∈ 0...2^{N-n}を求める
def seed_to_index0(s, n)
  if n == N
    0
  elsif ((s >> n) & 1) != 0
    (seed_to_index0(step_seed(s, 2**n), n + 1) * 2 - 1) % 2**(N-n)
  else
    seed_to_index0(s, n + 1) * 2
  end
end

前回の記事では下から5行目がseed_to_index0(step_seed(s, -2**n), n + 1) * 2 + 1となっていたのだが、-2^n進める定数をあらかじめ表にしておくより、2^n進める定数を表にしておく方が自然なので今回の方が好みだ。


トラックバック - http://d.hatena.ne.jp/oupo/20121122/1353554591