今さら森博嗣「笑わない数学者」のビリヤードの問題を解く(6)～答え合わせ編

0. はじめに

前回の記事までにビリヤードの問題は「ゴロム分度器」または「魔円陣」として古くから知られる数学の問題であることが分かった。今回，「魔円陣」を主題にした書籍（参考文献参照）を入手したので答え合わせを行いたい。

• 今さら森博嗣「笑わない数学者」のビリヤードの問題を解く(1)～導入編 - Qiita
• 今さら森博嗣「笑わない数学者」のビリヤードの問題を解く(2)～計算機を使う - Qiita
• 今さら森博嗣「笑わない数学者」のビリヤードの問題を解く(3)～計算量オーダーを考える - Qiita
• 今さら森博嗣「笑わない数学者」のビリヤードの問題を解く(4)～世界記録に挑む - Qiita
• 今さら森博嗣「笑わない数学者」のビリヤードの問題を解く(5)～未解決問題に挑む - Qiita

1. ボールの個数と解の個数の関係

参考文献の書籍にビリヤード問題の解の個数の公式がまんま書いてあった。

素数 $p$，自然数 $m$ として，ボールの個数 $n = 1 + p^m$ の場合の解の個数 $s$ は次のようになる。

s = \frac{\phi\left( p^{2m} + p^m + 1 \right)}{6m} \tag{1}

ここで $\phi(x)$ はオイラーの φ 関数といい，$x$ 以下で $x$ と互いに素な自然数の個数を返す。

自然数 $x$ が素数 $p_1, p_2, p_3, \cdots , p_q$ のべき乗数の積として以下のように定義する。

\begin{aligned}
x &= p_1^{k_1} \times p_2^{k_2} \times p_3^{k_3} \times \cdots \times p_q^{k_q} \\
&= \prod_{i = 1}^{q} p_i^{k_i}
\end{aligned} \tag{2}

このときオイラーの φ 関数 $\phi(x)$ は次のように定義される。

\begin{aligned}
\phi(x) &= x \times \left( 1 - \frac{1}{p_1} \right) \times \left( 1 - \frac{1}{p_2} \right) \times \left( 1 - \frac{1}{p_3} \right) \times \cdots \times \left( 1 - \frac{1}{p_q} \right) \\
&= x \prod_{i = 1}^{q} \left( 1 - \frac{1}{p_i} \right)
\end{aligned} \tag{3}

とくに $x$ が素数の場合は $\phi(x) = x - 1$ となる。

ボールの個数 $n = 1 + p^m$ に対し，ボールの数値の総和 $\sigma = n(n - 1) + 1$ とおくと，式$(1)$は，

s = \frac{\phi(n(n - 1) + 1)}{6m} = \frac{\phi(\sigma)}{6m} \tag{4}

とも書ける。これまでボールの個数 $n$ が増えると概ね解の個数 $s$ も増える傾向が見られたが，不規則に増減を繰り返してもいるようにも見えた。ボールの数値の総和 $\sigma$ が素数か否かによって解の個数が変化していたことが分かる。

オイラー φ 関数の JavaScript によるサンプルコードを以下に示す。

オイラー φ 関数

function euler_totient(n) {
	var	r = n;
	for(var i = 2; i * i <= n; i++) {
		if(n % i != 0) continue;
		while(n % i == 0) n /= i;
		r -= r / i;
	}
	if(n > 1) r -= r / n;
	return r;
}

この公式を使ってボールの個数 129 個以下における解の個数を計算してみると，見事に自作プログラムで求めた結果と一致した。ついでにボールの個数 988 個の場合の解の個数も計算してみると 162,526 個である。コチラも十万個近くあると予想していたのでまあ正解である。

なお厳密に言うなら式 $(1)$ は「有限体の射影平面」を用いて得た場合の公式であるとのこと。おそらく，トリボナッチ数列を用いて完全差集合から求める手法も「有限体の射影平面」に基づいた手法と数学的に等価なのだろう。

表１　ボールの個数と解の個数の関係

ボールの個数 $n$	素数 $p$	指数 $m$	総和 $\sigma$	解の個数 $s$
3	2	1	7	1
4	3	1	13	2
5	2	2	3×7	1
6	5	1	31	5
8	7	1	3×19	6
9	2	3	73	4
10	3	2	7×13	6
12	11	1	7×19	18
14	13	1	3×61	20
17	2	4	3×7×13	6
18	17	1	307	51
20	19	1	3×127	42
24	23	1	7×79	78
26	5	2	3×7×31	30
28	3	3	757	42
30	29	1	13×67	132
32	31	1	3×331	110
33	2	5	7×151	30
38	37	1	3×7×67	132
42	41	1	1,723	287
44	43	1	3×631	210
48	47	1	37×61	360
50	7	2	3×19×43	126
54	53	1	7×409	408
60	59	1	3,541	590
62	61	1	3×13×97	384
65	2	6	3×19×73	72
68	67	1	3×7²×31	420
72	71	1	5,113	852
74	73	1	3×1,801	600
80	79	1	3×7²×43	588
82	3	4	7×13×73	216
84	83	1	19×367	1,098
90	89	1	8,011	1,335
98	97	1	3×3,169	1,056
102	101	1	10,303	1,717
104	103	1	3×3,571	1,190
108	107	1	7×13×127	1,512
110	109	1	3×7×571	1,140
114	113	1	13×991	1,980
122	11	2	3×7×19×37	648
126	5	3	19×829	828
128	127	1	3×5,419	1,806
129	2	7	7²×337	336
988	987	1	975,157	162,526

ボールの数値の総和 $\sigma$ が素数の場合はそのまま，合成数の場合は素因数分解して示している。以前作成した素因数分解プログラムが役に立った。

2. 他の解の導き方

これもまた参考文献の書籍に書いてある。一つでも解を導くことができれば，他の解も簡単に導ける。これは全パラメータを網羅して探索する場合に有効であろう。最初の解を発見出来たらプログラムを終了し，他の解は本章の手法を用いて探したほうが効率的である。全探索する場合，ボールの個数 $n$ に対して計算量は $O(10^n)$ だったり，良くても $O(n^5)$ だったりするが，本手法は $O(n^3)$ 程度で済むからだ。

ボールの個数 $n = 4$ のときの解の一つ $b_i = \lbrace 1, 2, 6, 4 \rbrace$ で具体例を示そう。ボールの数値の総和 $\sigma = n(n - 1) + 1 = 13$ である。

まず，ボールの配置 $b_i$ を完全差集合 $c_i$ に変換する。完全差集合 $c_i$ はボールの配置 $b_i$ の積算値であり，$c_i = \lbrace 0, 1, 3, 9, 13 \rbrace$ となる。

次に完全差集合の各要素 $c_i$ を $m$（$1 \le m < \sigma$）倍して $\sigma$ で割った剰余を求め，各要素を昇順にソートする。

表２　ボールの個数 $n = 4$ のとき：$(m \times c_i) \bmod \sigma$

$m$	ソート前					$\Longrightarrow$	$m$	ソート後
1	0	1	3	9	13		1	0	1	3	9	13
2	0	2	6	5	13		2	0	2	5	6	13
3	0	3	9	1	13		3	0	1	3	9	13
4	0	4	12	10	13		4	0	4	10	12	13
5	0	5	2	6	13		5	0	2	5	6	13
6	0	6	5	2	13		6	0	2	5	6	13
7	0	7	8	11	13		7	0	7	8	11	13
8	0	8	11	7	13		8	0	7	8	11	13
9	0	9	1	3	13		9	0	9	1	3	13
10	0	10	4	12	13		10	0	4	10	12	13
11	0	11	7	8	13		11	0	7	8	11	13
12	0	12	10	4	13		12	0	4	10	12	13

ソート後の各要素に対して差分を取り，ボール配置に戻してやる。回転対称形と左右対称形を除くようボール配置を正規化すると，別解 $b_i = \lbrace 1, 3, 2, 7 \rbrace$ が現れる。

表３　ボールの個数 $n = 4$ のとき：ボール配置

$m$	正規化前				$\Longrightarrow$	$m$	正規化後
1	1	2	6	4		1	1	2	6	4
2	2	3	1	7		2	1	3	2	7
3	1	2	6	4		3	1	2	6	4
4	4	6	2	1		4	1	2	6	4
5	2	3	1	7		5	1	3	2	7
6	2	3	1	7		6	1	3	2	7
7	7	1	3	2		7	1	3	2	7
8	7	1	3	2		8	1	3	2	7
9	1	2	6	4		9	1	2	6	4
10	4	6	2	1		10	1	2	6	4
11	7	1	3	2		11	1	3	2	7
12	4	6	2	1		12	1	2	6	4

このケースでは全ての $m$ に対して重複解を含めて正しい解が得られたが，一般的にはそうとは限らない。例えばボールの個数 $n = 5$ のときの唯一解 $b_i = \lbrace 1, 3, 10, 2, 5 \rbrace$ について同様に行ってみる。ボールの数値の総和 $\sigma = n(n - 1) + 1 = 21$ となる。

完全差集合 $c_i = \lbrace 0, 1, 4, 14, 16, 21 \rbrace$ に対して，各要素 $c_i$ を $m$（$1 \le m < \sigma$）倍して $\sigma$ で割った剰余を求め，各要素を昇順にソートする。

表４　ボールの個数 $n = 5$ のとき：$(m \times c_i) \bmod \sigma$

$m$	ソート前						$\Longrightarrow$	$m$	ソート後
1	0	1	4	14	16	21		1	0	1	4	14	16	21
2	0	2	8	7	11	21		2	0	2	7	8	11	21
3	0	3	12	0	6	21		3	0	0	3	6	12	21
4	0	4	16	14	1	21		4	0	1	4	14	16	21
5	0	5	20	7	17	21		5	0	5	7	17	20	21
6	0	6	3	0	12	21		6	0	0	3	6	12	21
7	0	7	7	14	7	21		7	0	7	7	7	14	21
8	0	8	11	7	2	21		8	0	2	7	8	11	21
9	0	9	15	0	18	21		9	0	0	9	15	18	21
10	0	10	19	14	13	21		10	0	10	13	14	19	21
11	0	11	2	7	8	21		11	0	2	7	8	11	21
12	0	12	6	0	3	21		12	0	0	3	6	12	21
13	0	13	10	14	19	21		13	0	10	13	14	19	21
14	0	14	14	7	14	21		14	0	7	14	14	14	21
15	0	15	18	0	9	21		15	0	0	9	15	18	21
16	0	16	1	14	4	21		16	0	1	4	14	16	21
17	0	17	5	7	20	21		17	0	5	7	17	20	21
18	0	18	9	0	15	21		18	0	0	9	15	18	21
19	0	19	13	14	10	21		19	0	10	13	14	19	21
20	0	20	17	7	5	21		20	0	5	7	17	20	21

ソート後の各要素に対して差分を取り，ボール配置に戻してやる。ただし，ボール配置を正規化するため１の玉を探そうとしても見つからないNGケースが８個ある。NGにならない，すなわち重複解も含めて解が得られるケースは１２個あるが，おそらくこれは $\phi(21) = 12$ だからだと考えられる。ちなみに先のケースも $\phi(13) = 12$ だったので全てのケースで解が得られたのだろう。

表５　ボールの個数 $n = 5$ のとき：ボール配置

$m$	正規化前					$\Longrightarrow$	$m$	正規化後
1	1	3	10	2	5		1	1	3	10	2	5
2	2	5	1	3	10		2	1	3	10	2	5
3	0	3	3	6	9		3	NG
4	1	3	10	2	5		4	1	3	10	2	5
5	5	2	10	3	1		5	1	3	10	2	5
6	0	3	3	6	9		6	NG
7	7	0	0	7	7		7	NG
8	2	5	1	3	10		8	1	3	10	2	5
9	0	9	6	3	3		9	NG
10	10	3	1	5	2		10	1	3	10	2	5
11	2	5	1	3	10		11	1	3	10	2	5
12	0	3	3	6	9		12	NG
13	10	3	1	5	2		13	1	3	10	2	5
14	7	7	0	0	7		14	NG
15	0	9	6	3	3		15	NG
16	1	3	10	2	5		16	1	3	10	2	5
17	5	2	10	3	1		17	1	3	10	2	5
18	0	9	6	3	3		18	NG
19	10	3	1	5	2		19	1	3	10	2	5
20	5	2	10	3	1		20	1	3	10	2	5

どんな $m$ をかけても解が得られるという訳ではないことが分かったと思う。正しくは「有限体の生成元」をかければ良いらしいが，その値を求めるよりも「$m = 2, 3, 4, 5, \cdots $ と順に試したほうが案外早く見つかる」と参考文献にも書いてある。

ボールの数値の総和 $\sigma$ が素数の場合，すべての $m$ が生成元になり得るが，合成数の場合，約数とその倍数は生成元になり得ないような感じである。たとえば $\sigma = 21$ のとき，$m$ が３の倍数と７の倍数のときは解を得られていない。

実装コードを以下に示す。今回も JavaScript（WSH）で組んでみた。

AnotherSolver.JS

var	args = WScript.Arguments.Unnamed;
var	ret = main(args)
try {
	WScript.Quit(ret);
} catch(e) {
	/* 何もしない */
}
//------------------------------------------------------------------------------
// メイン関数
//------------------------------------------------------------------------------
function main(args) {
	//--------------------------------------------------------------------------
	// ヘルプメッセージ
	//--------------------------------------------------------------------------
	if(args.Count == 0) {
		WScript.Echo("Usage: AnotherSolver(.js) (option) [ball-position] ...");
		return false;
	}
	//--------------------------------------------------------------------------
	// コマンドライン解析
	//--------------------------------------------------------------------------
	var	ball = [];
	for(var i = 0; i < args.Count; i++) {
		var	n = parseInt(args(i));
		ball.push(n);
	}
	//--------------------------------------------------------------------------
	// ボール配置のチェック
	//--------------------------------------------------------------------------
	if(!check(ball)) {
		WScript.Echo("ボール配置が不正です!!");
		return false;
	}
	//--------------------------------------------------------------------------
	// ボール配置を完全差集合に変換する
	//--------------------------------------------------------------------------
	var	pcds = [];
	for(var i = 0, n = 0; i < ball.length; i++) {
		pcds.push(n);
		n += ball[i];
	}
	var	ball_sum = n;
	//--------------------------------------------------------------------------
	// 別解を探す
	//--------------------------------------------------------------------------
	var	d = [];
	for(var m = 1; m < ball_sum; m++) {
		//----------------------------------------------------------------------
		// 差集合を生成
		//----------------------------------------------------------------------
		var	s = [];
		for(var i = 0; i < pcds.length; i++) {
			var	n = pcds[i] * m % ball_sum; 
			s.push(n);
		}
		s.push(ball_sum);
		s.sort(function(a, b) {
			return a - b;
		});
		//----------------------------------------------------------------------
		// ボール配置に変換
		//----------------------------------------------------------------------
		var	t = [];
		for(var i = 1; i < s.length; i++)
			t.push(s[i] - s[i - 1]);
		//----------------------------------------------------------------------
		// 回転対称形を除く（１のボールを頭出し）
		//----------------------------------------------------------------------
		for(var i = 0; i < t.length; i++)
			if(t[i] == 1) break;
		if(i >= t.length) continue;
		var	u = [];
		for(var j = i, k = 0; k < t.length; k++) {
			u.push(t[j]);
			if(++j >= t.length) j = 0;
		}
		//----------------------------------------------------------------------
		// 左右対称形を除く
		//----------------------------------------------------------------------
		u.push(1);
		if(u[1] > u[u.length -2]) u.reverse();
		u.pop();
		//----------------------------------------------------------------------
		// ボール配置のチェック＆重複するボール配置を除く
		//----------------------------------------------------------------------
		if(check(u)) {
			var	v = u.join(" ");
			if(d[v] == null) d[v] = m;
		}
	}
	//--------------------------------------------------------------------------
	// ボール配置をソートして出力する
	//--------------------------------------------------------------------------
	var	e = [];
	for(var v in d)
		e.push({ key:d[v], val:v, arr:v.split(" ") });
	e.sort(function(a, b) {
		for(var i = 0; i < a.arr.length; i++)
			if(a.arr[i] != b.arr[i]) return a.arr[i] - b.arr[i];
		return 0;
	});
	for(var i = 0; i < e.length; i++)
		WScript.Echo([e[i].key, e[i].val].join("\t"));
	//--------------------------------------------------------------------------
	// 解の個数を出力する
	//--------------------------------------------------------------------------
	WScript.Echo(e.length);
	return true;
}
//------------------------------------------------------------------------------
// ボール配置のチェック
//------------------------------------------------------------------------------
function check(ball) {
	var	num = ball.length;
	var	sum = num * (num -  1) + 1;
	for(var i = 0, s = 0; i < ball.length; i++)
		s += ball[i];
	if(sum != s) return false;
	var	flag = [];
	for(var i = 0; i <= sum; i++)
		flag.push(false);
	for(var i = 0; i < ball.length - 1; i++) {
		for(var j = 0, s = 0; i + j >= 0; j--) {
			s += ball[i + j];
			if(flag[s] || flag[sum - s]) return false;
			flag[s] = true;
			flag[sum - s] = true;
		}
	}
	return true;
}

実行例を以下に示す。ボールの個数 $n = 17$ のときの解の一つを与えると，その解を含めたすべての解を求めることができる。各行一つ目の数値は倍数の $m$ であり，続く数列がボールの配置である。最後の行の 6 は解の個数である。

実行例

anothersolver 1 2 4 8 16 32 27 26 11 9 45 13 10 29 5 17 18
1       1 2 4 8 16 32 27 26 11 9 45 13 10 29 5 17 18
23      1 3 12 10 31 7 27 2 6 5 19 20 62 14 9 28 17
29      1 7 3 15 33 5 24 68 2 14 6 17 4 9 19 12 34
19      1 7 12 44 25 41 9 17 4 6 22 33 13 2 3 11 23
5       1 7 31 2 11 3 9 36 17 4 22 6 18 72 5 10 19
11      1 21 11 50 39 13 6 4 14 16 25 26 3 2 7 8 27
6

上記のスクリプトはボールの個数 100 個前後であれば問題ないが，1000 個くらいになると内部メモリがオーバーフローして異常終了してしまうので注意されたい。

3. まとめ

最後に残された謎について書こうと思ったのだが，本記事は長くなり過ぎたので次回にする。

4. 参考文献

魔円陣を主題にした書籍は少ないので貴重である。

オイラーの φ 関数の定義と解説はコチラ