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

• 前回の記事: 線形合同法であるseedが0からいくつ進めたものかを得る - oupoの日記
• 続きの記事: 線形合同法であるseedが0からいくつ進めたものかを得る その3 (法が2のべきとは限らない一般の場合) - oupoの日記

周期2^k, 初項0のLCG数列{x_n}に対して、{x_{2n} / 2}って周期2^(k-1)のLCG数列だよね、というアイディアで再実装してみました。

#include <stdio.h>
#include <stdint.h>

typedef uint32_t u32;

/* x[0] = 0, x[n+1] = (a * x[n] + b) mod 2^k
   で与えられる数列{x[n]}において、x[n] = sを満たすnを求める
   ただし数列の周期は2^kになっていなければならない */
u32 calc_index(u32 a, u32 b, u32 s, u32 k) {
	if (k == 0) {
		return 0;
	} else if ((s & 1) == 0) {
		return calc_index(a * a, (a + 1) * b / 2, s / 2, k - 1) * 2;
	} else {
		return calc_index(a * a, (a + 1) * b / 2, (a * s + b) / 2, k - 1) * 2 - 1;
	}
}

int main(void) {
	u32 A = 0x41c64e6d;
	u32 B = 0x6073;

	printf("%.8x\n", calc_index(A, B, 0, 32));
	printf("%.8x\n", calc_index(A, B, B, 32));
	printf("%.8x\n", calc_index(A, B, A * B + B, 32));
	printf("%.8x\n", calc_index(A, B, A * A * B + A * B + B, 32));
	return 0;
}