数列
xiが例えば32ビットの語の範囲の値しかとれなくて,
xi+1=
f(
xi)を繰り返していくと, いつかは昔のある値が現れ, 以後同じ数列を辿ることになって, ループを回り始める.
このループの開始位置とループ長を求めるのがループ検出である.
そのひとつの方法にBill Gosperの考案したものがあり, あちこちに説明があるが, あまり真面目に読んだことはなかった. しかし気にはなっていたので, この数日プログラムを書いてみたりしながら, その仕掛けを調べてみた.
いま
x1,
x2,
x3, ... が1,2,3,4,5,6,7,8,9.10,6,7,8,9,10,6,7,8,...であったとしよう.(
x0から始める流儀もあり, また私はそれが好きだが, 今回の添字は1から始めることにした.)
The Art of Computer Programming(TAOCP)の演習問題3.1,6の記述を借用すれば, 「
x1,
x2,...,
xμ,...,
xμ+λ-1は互いに相異なるが,
n ≥μの時は
xn=
xn+λであるμとλがある.」
このμが開始位置, λがループ長で, 上の例ではμ=6, λ=5だ.
記憶場所も時間もふんだんにある時は,
xの値をすべて順に記憶しておき, 新しい
xについて, それを過去のすべての
xの値のそれぞれと比較する. そこに同じものがあれば, 前の値のあったところがループの開始場所であり, そこから今の値までの長さがループ長である.
しかし, こんなことはとてもやっていられないから, あの手この手を考えることになる. その一つにGosperの提案したものがある. 以下は添字などを私流に多少修正したものである.
これは簡単にいえば, 繰り返しが「現れるや否や」発見するのは諦める, 最新の
xの1つ前, 2つ前, 4つ前, 8つ前, ... などを覚えておき, ループに入り同じ値が現れ始めたらなるべく早目にそれを発見しようという方針である.
「1つ, 2つ, 4つ, 8つ,...」から推察されるように,
n個の
xを調べた時は, log
2 n個程度の配列
Tを使う.
ところで,
n>0を二進法で表した時, もっとも右の1の位置, もっとも左の1の位置をTAOCPに倣い, それぞれρ(
n), &lambda(
n)とする. ビットの位置は右端を0とする. (上のループ長の定義にもλがあって申し訳ない.)
ρ(1)=0, λ(1)=0;
ρ(2)=1, λ(2)=1;
ρ(3)=0, λ(3)=1;
ρ(4)=2, λ(4)=2;
ρ(5)=0, λ(5)=2;
ρ(6)=1, λ(6)=2;
ρ(7)=0, λ(7)=2;
ρ(8)=3, λ(8)=3;
ρの方は別名ルーラー関数(物差関数)という.
十進法の物差しとは様子が違い, ひと目盛りが半分ずつに分割されていて, インチの物差しが大体こういう感じである.
ρとλを計算する関数 rhoとlamは次のようだ.
ρは
(define (rho x)
(define (r x rh)
(if (= (modulo x 2) 0) (r (quotient x 2) (+ rh 1)) rh))
(if (= x 0) 'error (r x 0)))
(map rho (a2b 1 16)) => (0 1 0 2 0 1 0 3 0 1 0 2 0 1 0)
またλは関数logがあれば
(define (lam n)
(inexact->exact (floor (/ (log n) (log 2))))
(map lam (a2b 1 16))=>
(0 1 1 2 2 2 2 3 3 3 3 3 3 3 3)
ρのように自力で計算するには
(define (lam x)
(cond ((= x 0) 'error)
((= x 1) 0)
(else (+ (lam (quotient x 2)) 1))))
さて配列をT0, T1, T2, ... とする. (こちらは0から始める) そしてxiをTρ(i)に順に置くのである.
xiが次々と置かれていく様子は次のようだ.
| i | T0 |
T1 |
T2 |
T3 |
| 1 | x1 |
| 2 | x1 |
x2 |
| 3 | x3 |
x2 |
| 4 | x3 |
x2 | x4 |
| 5 | x5 |
x2 | x4 |
| 6 | x5 |
x6 | x4 |
| 7 | x7 |
x6 | x4 |
| 8 | x7 |
x6 | x4 |
x8 |
| 9 | x9 |
x6 | x4 |
x8 |
| 10 | x7 |
x10 | x4 |
x8 |
注意すべきは, この表のうち見えているのは最後の行だけということである. つまりi=8の時には, x1, | x2, | x3, | x5はすでに消えている.
さて, 上の例のようなループがあった時の同様な配列のx10=10までの変化の様子は,
| i | T0 |
T1 |
T2 |
T3 |
| 1 | 1 |
| 2 | 1 |
2 |
| 3 | 3 |
2 |
| 4 | 3 |
2 | 4 |
| 5 | 5 |
2 | 4 |
| 6 | 5 |
6 | 4 |
| 7 | 7 |
6 | 4 |
| 8 | 7 |
6 | 4 |
8 |
| 9 | 9 |
6 | 4 |
8 |
| 10 | 9 |
10 | 4 |
8 |
のようだ. ここでx11=6になった時, この配列は9, 10, 4, 8であり, この中に6はないから, 6を(ρ(11)=0だから)T0のところに置く. したがって配列は6, 10, 4, 8になる. 次, x12=7も配列にないから, 7をT2に置き, 配列は6, 10, 7, 8になる.
次はx13=8だが, i=8の時に作ったT3に8があるので, 循環が始まっていたのがわかる.
そこで問題はiの行のTkを置いたのは誰かということだ. 最初の表からxを削除した表は次のようだが
| i\k | 0 |
1 |
2 |
3 |
| 1 | 1 |
| 2 | 1 |
2 |
| 3 | 3 |
2 |
| 4 | 3 |
2 | 4 |
| 5 | 5 |
2 | 4 |
| 6 | 5 |
6 | 4 |
| 7 | 7 |
6 | 4 |
| 8 | 7 |
6 | 4 |
8 |
| 9 | 9 |
6 | 4 |
8 |
| 10 | 7 |
10 | 4 |
8 |
これで例えばi=7, k=2ならこの4を置いたのはi=4であったということが知りたい.
いろいろやってみると, iからiの下からk+1ビットをとったものの, 先頭の0,1を逆にしたものを引くといいらしいことが分かる.
iの下からk+1ビットをとるには, 2k+2-1で1のビット列を作ってマスクする.
i &(2k+2-1)
これを 2kでxorすると先頭のビットが反転する.
それをiから引くのである.
Schemeで書くと
(define (m i k)
(- i (fix:xor (fix:and i (- (fix:lsh 1 (+ k 1)) 1))
(fix:lsh 1 k))))
であって, 実際に計算すると
(do ((i 1 (+ i 1))) ((= i 11))
(display (cons i
(map (lambda (k) (m i k)) (a2b 0 (+ (lam i) 1)))))
(newline))
(1 1)
(2 1 2)
(3 3 2)
(4 3 2 4)
(5 5 2 4)
(6 5 6 4)
(7 7 6 4)
(8 7 6 4 8)
(9 9 6 4 8)
(10 9 10 4 8)
となって上の表と一致する.
従ってi=13の時のx13>=8がi=12, k=3にあり, (m 12 3)=8が分ったから, λ=13-8=5であった.
プログラムを書いてやってみよう. 最初の数列xiは, μとλから
(define (f i mu lm)
(if (< i mu) i (+ (modulo (- i mu) lm) mu)))
(map (lambda (i) (f i 6 5)) (a2b 1 19)) =>
(1 2 3 4 5 6 7 8 9 10 6 7 8 9 10 6 7 8)
と作る. プログラム全体と実行結果を下に示す.
(define (iloop i x)
(define (kloop k)
(if (> k (lam (- i 1))) (iloop i x)
(begin
(if (= x (vector-ref tab k)) (list i k (m (- i 1) k))
(kloop (+ k 1))))))
(vector-set! tab (rho i) x)
(display i) (display (take (+ (lam i) 1) (vector->list tab)))
(newline)
(set! i (+ i 1))
(set! x (f i mu lm))
(kloop 0))
(define mu 6) (define lm 5)
(iloop 1 (f 1 mu lm))
1(1)
2(1 2)
3(3 2)
4(3 2 4)
5(5 2 4)
6(5 6 4)
7(7 6 4)
8(7 6 4 8)
9(9 6 4 8)
10(9 10 4 8)
11(6 10 4 8)
12(6 10 7 8)
=> (13 3 8)
つまりi=13のとき, k=3の所に同じものがあり, その時のiは8であった. 従ってλは13-8=5であった.
| 
