LoginSignup
2
1

今さら森博嗣「笑わない数学者」のビリヤードの問題を解く(2)~計算機を使う

Last updated at Posted at 2024-01-28

0. 前回までのあらすじ

森博嗣氏の小説「笑わない数学者」に出てくる通称「ビリヤードの問題」は,人間が手計算でも解け,かつ解が唯一定まる(対称形を除く)美しい問題であった。詳しくは前回の記事を参照して欲しい。

1. 一般化した問題

$n$ 個のボールを円周上に並べる。ボールには自然数が書かれており,隣り合うボールを取り出し,取り出したボールの数値の合計が $1$ から $\sigma = n(n - 1) + 1$ までの全ての数を作ることができる場合のボールの配置を求めよ。

円周上に並べた $n$ 個のボールから隣り合うボールを取り出す組み合わせは,$1$~$n - 1$ 個を取り出す組み合わせが $n$ 通りずつ,$n$ 個全部を取り出す組み合わせが $1$ 通りなので,合計 $\sigma = n(n - 1) + 1$ 通りある。取り出したボールの数値の和が $1$ から $\sigma$ までの全ての数値を作れるということは,取り出したボールの数値が重複しないことを意味する。

2. インクリメンタル重複チェック

すべてのボールを配置し終えてから重複チェックを行うと時間がかかり,エラー時の損失も大きいことからインクリメンタルな重複チェック,すなわちボールを一つ配置する毎に可能な範囲での重複チェックを行うことを考える。$k$ 番目のボール $b_k$ を配置するとき,$b_k$ も含めて以前に配置したボールの数値 $b_i$($1 \le i \le k$)しか得られないが,総和が $\sigma$ になることを利用すれば早期の重複チェックを行うことができる。以下は $n = 5$ の場合の例である。前回の記事で述べたように $n = 5$ の場合は $21$ 通りの組み合わせがある。チェックの順番自体はシャッフルされているが,全ての組み合わせを網羅できている。

表1 $n = 5$ の場合
$b_1$ $b_2$ $b_3$ $b_4$ $b_5$ 数値 チェック時期
        $b_1$ $b_1$ 配置時
  $\sigma - b_1$
        $b_2$ $b_2$ 配置時
  $\sigma - b_2$
      $b_2 + b_1$
    $\sigma - (b_2 + b_1)$
        $b_3$ $b_3$ 配置時
  $\sigma - b_3$
      $b_3 + b_2$
    $\sigma - (b_3 + b_2)$
    $b_3 + b_2 + b_1$
      $\sigma - (b_3 + b_2 + b_1)$
        $b_4$ $b_4$ 配置時
  $\sigma - b_4$
      $b_4 + b_3$
  $\sigma - (b_4 + b_3)$
    $b_4 + b_3 + b_2$
      $\sigma - (b_4 + b_3 + b_2)$
  $b_4 + b_3 + b_2 + b_1$
        $\sigma - (b_4 + b_3 + b_2 + b_1)$
$\sigma$ 不要

3. 探索範囲を絞る

3.1 1番目のボール

回転対称形を除外して,$1$ 番目に置くボールの数値は $b_1 = 1$ とする。

3.2 2番目のボールの範囲

左右対称形を除外して,$2$ 番目に置くボールの数値 $b_2$ は $b_1$ の逆隣りの $b_n$ より小さい,すなわち $b_2 < b_n$ とする。

\sigma = 1 + b_2 + \underbrace{\sum_{i = 3}^{n - 1} b_i}_{(*1)} + b_n

$(*1)$ の最小値は

(*1) = \sum_{i = 3}^{n - 1} b_i \ge \sum_{k = 2}^{n - 2} k = \frac{n(n - 3)}{2}

であり,$b_n$ の最小値は

b_n \ge b_2 + 1

であるから

\sigma = 1 + b_2 + \sum_{i = 3}^{n - 1} b_i + b_n \ge 2 + 2b_2 + \frac{n(n - 3)}{2}

であり

b_2 \le \left\lfloor \frac{(n + 2)(n - 1)}{4} \right\rfloor

という条件が導ける。

3.3 3番目以降のボールの範囲

$3$ 番目以降に置くボールの数値 $b_k$($3 \le k < n$)の範囲について考える。

\sigma = \underbrace{\sum_{i = 1}^{k - 1} b_i}_{(*1)} + b_k + \underbrace{\sum_{i = k + 1}^{n - 1}  b_i}_{(*2)} + b_n

ここで $(*2)$ の最小値は $(*1)$ で使用していない $2$ 以上の自然数の和であり,$b_n$ の最小値は $(*1)$ および $(*2)$ で使用していない $b_2$ より大きい自然数となる。

4. コード解説

この問題は $n$ が大きくなるとガチで計算機パワーを必要とするので,ガチガチに最適化のできる C 言語で組む。

4.1 変数定義

C 言語の配列は 0-origin なので,分かり易く 1-origin なるようにした。すなわち $1$ 番目のボールの数値は ball[1] である。重複フラグテーブル flag[s]true のとき,取り出したボールの数値の総和 s が既に使用中であることを意味し,false のときは未使用とする。

変数宣言部
static	int		count;						// 解の数
static	int		ball_num;					// ボールの個数
static	int		ball_sum;					// ボールの総和
static	char	ball[1 + MAX_BALL_NUM];		// ボールの配置テーブル
static	bool	flag[1 + MAX_BALL_SUM];		// 重複フラグテーブル true:使用,false:未使用

4.2 ボールの範囲を求める

$2$ 番目のボールの範囲は簡単に求まる。

2番目のボールの範囲を求める
range = ( ball_num + 2 ) * ( ball_num - 1 ) / 4;

$3$ 番目以降のボールの範囲は簡単には求まらない。sumball[1] から ball[pos - 1] までの総和とする。sball[pos + 1] から ball[ball_num] までの総和の最小値とする。すなわち $2$ 以上の未使用の自然数の総和であり,かつ ball[ball_num]ball[2] よりも大きいことが条件となる。

3番目以降のボールの範囲を求める
i = 2, s = 0;
for( j = pos + 1; j < ball_num; j++ ) {
	while( flag[i] ) i++;
	s += i++;
}
if( i < ( j = ball[2] + 1 ) ) i = j;
while( flag[i] ) i++;
s += i;
range = ball_sum - sum - s;

4.3 重複チェック

ボールの数値を i とする。まず重複フラグテーブルを探索して未使用の数値を探す。

重複フラグのスキップ
while( flag[i] ) i++;
pos = i;

次に事前チェックを行う。重複フラグを書き換えてしまうとエラー時のペナルティ(フラグをリセットする処理)がそれなりに重いので,重複フラグを書き換えない簡易チェックを行って次のインクリメンタル重複チェックに進む事例を減らしている。直前のスキップ処理の結果より flag[pos]false であることから初回の判定を省略しており,また pos >= 3 であることから最初の三回分をループ展開している。ここは高頻度で実行される部分なので,小さな最適化を積み重ねている。

重複フラグを書き換えない事前チェック
j = pos; s = ball[j--];
s += ball[j--]; if( flag[s] ) goto SKIP;
s += ball[j--]; if( flag[s] ) goto SKIP;
while( j > 0 ) {
	s += ball[j--]; if( flag[s] ) goto SKIP;
}

その次にインクリメンタルな重複チェックを行う。エラー時には必要最小限,すなわちセットした分だけリセットを行う。ここはコードが複雑なのでループ展開しないほうが良いようだ。

インクリメンタル重複チェック
j = pos; s = ball[j];
for(;;) {
	if( flag[s] || flag[ball_sum - s] ) {
		while( j < pos ) {
			s -= ball[j++];
			flag[s] = false; flag[ball_sum - s] = false;
		}
		goto SKIP;
	}
	flag[s] = true; flag[ball_sum - s] = true;
	if( --j == 0 ) break;
	s += ball[j];
}

5. 実装コード

実装コードを以下に示す。ブール型を使っているので C99 に対応した C コンパイラならば通るだろう。自分は Visual Studio 2022 Community Edition を用いた。

billiard.c
#include <stdio.h>
#include <stdlib.h>
#include <stdbool.h>
//------------------------------------------------------------------------------
// 定数定義
//------------------------------------------------------------------------------
#define	MIN_BALL_NUM	(4)									// ボールの個数の最小値
#define	MAX_BALL_NUM	(20)								// ボールの個数の最大値
#define	MAX_BALL_SUM	(1+MAX_BALL_NUM*(MAX_BALL_NUM-1))	// ボールの総和の最大値
//------------------------------------------------------------------------------
// グローバル変数
//------------------------------------------------------------------------------
static	int		count;					// 解の数
static	int		ball_num;				// ボールの個数
static	int		ball_sum;				// ボールの総和
static	char	ball[1 + MAX_BALL_NUM];	// ボールの配置テーブル
static	bool	flag[1 + MAX_BALL_SUM];	// 重複フラグテーブル true:使用,false:未使用
//------------------------------------------------------------------------------
// 3rd ボール以降のチェック
// pos: 配置するボールの位置,3 ≦ pos < ball_num
// sum: ball[1] ~ ball[pos - 1] までの総和
//------------------------------------------------------------------------------
static	void	check_collision( int pos, int sum ) {
	int	i, j, s, range;
	//--------------------------------------------------------------------------
	// ボールの範囲を求める
	//--------------------------------------------------------------------------
	i = 2, s = 0;
	for( j = pos + 1; j < ball_num; j++ ) {
		while( flag[i] ) i++;
		s += i++;
	}
	if( i < ( j = ball[2] + 1 ) ) i = j;
	while( flag[i] ) i++;
	s += i;
	range = ball_sum - sum - s;
	//--------------------------------------------------------------------------
	// ボールの探索
	//--------------------------------------------------------------------------
	i = 2; while( flag[i] ) i++;
	while( i <= range ) {
		//----------------------------------------------------------------------
		// ボールの配置
		//----------------------------------------------------------------------
		ball[pos] = (char)i;
		//----------------------------------------------------------------------
		// 重複フラグを書き換えない事前チェック
		// pos ≧ 3 であることを利用してループ展開
		//----------------------------------------------------------------------
		j = pos; s = ball[j--];
		s += ball[j--]; if( flag[s] ) goto SKIP;
		s += ball[j--]; if( flag[s] ) goto SKIP;
		while( j > 0 ) {
			s += ball[j--]; if( flag[s] ) goto SKIP;
		}
		//----------------------------------------------------------------------
		// インクリメンタル重複チェック
		//----------------------------------------------------------------------
		j = pos; s = ball[j];
		for(;;) {
			if( flag[s] || flag[ball_sum - s] ) {
				//--------------------------------------------------------------
				// 重複フラグの部分リセット
				//--------------------------------------------------------------
				while( j < pos ) {
					s -= ball[j++];
					flag[s] = false; flag[ball_sum - s] = false;
				}
				goto SKIP;
			}
			//------------------------------------------------------------------
			// 重複フラグのセット
			//------------------------------------------------------------------
			flag[s] = true; flag[ball_sum - s] = true;
			if( --j == 0 ) break;
			s += ball[j];
		}
		//----------------------------------------------------------------------
		// 重複チェックをパスした場合
		//----------------------------------------------------------------------
		if( pos + 1 < ball_num ) {
			//------------------------------------------------------------------
			// 次のボールの重複チェック
			//------------------------------------------------------------------
			check_collision( pos + 1, sum + i );
		} else {
			//------------------------------------------------------------------
			// 解発見時,ボールの配置を出力する
			//------------------------------------------------------------------
			ball[pos + 1] = (char)( ball_sum - sum - i );
			fprintf( stdout, "(%d)", ++count );
			for( int k = 1; k <= ball_num; k++ )
				fprintf( stdout, " %d", ball[k] );
			fprintf( stdout, "\n" );
		}
		//----------------------------------------------------------------------
		// 重複フラグのリセット
		//----------------------------------------------------------------------
		j = pos; s = ball[j];
		for(;;) {
			flag[s] = false; flag[ball_sum - s] = false;
			if( --j == 0 ) break;
			s += ball[j];
		}
SKIP:
		i++; while( flag[i] ) i++;
	}
}
//------------------------------------------------------------------------------
// メイン関数
//------------------------------------------------------------------------------
int		main( int argc, char *argv[] ) {
	//--------------------------------------------------------------------------
	// ヘルプメッセージ
	//--------------------------------------------------------------------------
	if( argc < 2 ) {
		fprintf( stderr, "Usage: billiard(.exe) [ball-num]\n" );
		return( -1 );
	}
	//--------------------------------------------------------------------------
	// ボールの個数のチェック
	//--------------------------------------------------------------------------
	ball_num = atoi( argv[1] );
	if( ball_num < MIN_BALL_NUM || ball_num > MAX_BALL_NUM ) {
		fprintf( stderr, "ボールの個数 %d が不正です!!\n", ball_num );
		return( -1 );
	}
	//--------------------------------------------------------------------------
	// ボールの数値の総和を求める
	//--------------------------------------------------------------------------
	ball_sum = 1 + ball_num * ( ball_num - 1 );
	//--------------------------------------------------------------------------
	// 1st ボールの配置
	//--------------------------------------------------------------------------
	ball[1] = 1;
	//--------------------------------------------------------------------------
	// 重複フラグのセット
	//--------------------------------------------------------------------------
	flag[1] = true; flag[ball_sum - 1] = true;
	//--------------------------------------------------------------------------
	// 2nd ボールの範囲を求める
	//--------------------------------------------------------------------------
	int	range = ( ball_num + 2 ) * ( ball_num - 1 ) / 4;
	//--------------------------------------------------------------------------
	// 2nd ボールの探索
	//--------------------------------------------------------------------------
	for( int i = 2; i <= range; i++ ) {
		//----------------------------------------------------------------------
		// 2nd ボールの配置
		//----------------------------------------------------------------------
		ball[2] = (char)i;
		//----------------------------------------------------------------------
		// 重複フラグのセット
		//----------------------------------------------------------------------
		flag[i]     = true; flag[ball_sum - i]     = true;
		flag[i + 1] = true; flag[ball_sum - i - 1] = true;
		//----------------------------------------------------------------------
		// 3rd ボール以降の重複チェック
		//----------------------------------------------------------------------
		check_collision( 3, i + 1 );
		//----------------------------------------------------------------------
		// 重複フラグのリセット
		//----------------------------------------------------------------------
		flag[i]     = false; flag[ball_sum - i]     = false;
		flag[i + 1] = false; flag[ball_sum - i - 1] = false;
	}
	return( 0 );
}

C 言語のシンタクスハイライトはブール型にまだ対応してないようなので(2024年1月時点)C++ を指定してみた。

6. 実行例

本プログラムはコマンドプロンプトから実行する。引数としてボールの個数を指定する。以下は $n = 12$ の場合である。最近のプロセッサであれば実行時間は十秒かからないと思う。

実行例
billiard 12
(1) 1 2 9 8 14 4 43 7 6 10 5 24
(2) 1 2 12 31 25 4 9 10 7 11 16 5
(3) 1 2 14 4 37 7 8 27 5 6 13 9
(4) 1 2 14 12 32 19 6 5 4 18 13 7
(5) 1 3 8 9 5 19 23 16 13 2 28 6
(6) 1 3 12 34 21 2 8 9 5 6 7 25
(7) 1 3 23 24 6 22 10 11 18 2 5 8
(8) 1 4 7 3 16 2 6 17 20 9 13 35
(9) 1 4 16 3 15 10 12 14 17 33 2 6
(10) 1 4 19 20 27 3 6 25 7 8 2 11
(11) 1 4 20 3 40 10 9 2 15 16 6 7
(12) 1 5 12 21 29 11 3 16 4 22 2 7
(13) 1 7 13 12 3 11 5 18 4 2 48 9
(14) 1 8 10 5 7 21 4 2 11 3 26 35
(15) 1 14 3 2 4 7 21 8 25 10 12 26
(16) 1 14 10 20 7 6 3 2 17 4 8 41
(17) 1 15 5 3 25 2 7 4 6 12 14 39
(18) 1 22 14 20 5 13 8 3 4 2 10 31

7. 計算結果まとめ

$n \le 17$ までの計算結果を以下に示す。$n$ が大きくなるにつれて,解の数が爆発的に増えていくものと思われたが,意外とそうでもなさそうだ。むしろ $n$ が大きくなると解が存在しない例が多くなるような気もする。残念ながら自分の数学の能力では解が存在する条件は導けなかった。

表2 計算結果まとめ(計算時間の単位は秒)
$n$ 計算時間 解の数 ボール配置
4 0.01 2 1 2 6 4
1 3 2 7
5 0.01 1 1 3 10 2 5
6 0.01 5 1 2 5 4 6 13
1 2 7 4 12 5
1 3 2 7 8 10
1 3 6 2 5 14
1 7 3 2 4 14
7 0.01 0 なし
8 0.01 6 1 2 10 19 4 7 9 5
1 3 5 11 2 12 17 6
1 3 8 2 16 7 15 5
1 4 2 10 18 3 11 8
1 4 22 7 3 6 2 12
1 6 12 4 21 3 2 8
9 0.03 4 1 2 4 8 16 5 18 9 10
1 4 7 6 3 28 2 8 14
1 6 4 24 13 3 2 12 8
1 11 8 6 4 3 2 22 16
10 0.15 6 1 2 6 18 22 7 5 16 4 10
1 3 9 11 6 8 2 5 28 18
1 4 2 20 8 9 23 10 3 11
1 4 3 10 2 9 14 16 6 26
1 5 4 13 3 8 7 12 2 36
1 6 9 11 29 4 8 2 3 18
11 1.17 0 なし
12 11.39 18 1 2 9 8 14 4 43 7 6 10 5 24
1 2 12 31 25 4 9 10 7 11 16 5
1 2 14 4 37 7 8 27 5 6 13 9
1 2 14 12 32 19 6 5 4 18 13 7
1 3 8 9 5 19 23 16 13 2 28 6
1 3 12 34 21 2 8 9 5 6 7 25
1 3 23 24 6 22 10 11 18 2 5 8
1 4 7 3 16 2 6 17 20 9 13 35
1 4 16 3 15 10 12 14 17 33 2 6
1 4 19 20 27 3 6 25 7 8 2 11
1 4 20 3 40 10 9 2 15 16 6 7
1 5 12 21 29 11 3 16 4 22 2 7
1 7 13 12 3 11 5 18 4 2 48 9
1 8 10 5 7 21 4 2 11 3 26 35
1 14 3 2 4 7 21 8 25 10 12 26
1 14 10 20 7 6 3 2 17 4 8 41
1 15 5 3 25 2 7 4 6 12 14 39
1 22 14 20 5 13 8 3 4 2 10 31
13 112.67 0 なし
14 1124.96 20 1 2 13 7 5 14 34 6 4 33 18 17 21 8
1 2 21 17 11 5 9 4 26 6 47 15 12 7
1 2 28 14 5 6 9 12 48 18 4 13 16 7
1 3 5 6 25 32 23 10 18 2 17 7 22 12
1 3 12 7 20 14 44 6 5 24 2 28 8 9
1 3 24 6 12 14 11 55 7 2 8 5 16 19
1 4 6 31 3 13 2 7 14 12 17 46 8 19
1 4 8 52 3 25 18 2 9 24 6 10 7 14
1 4 20 2 12 3 6 7 33 11 8 10 35 31
1 5 2 24 15 29 14 21 13 4 33 3 9 10
1 5 23 27 42 3 4 11 2 19 12 10 16 8
1 6 8 22 4 5 33 21 3 20 32 16 2 10
1 8 3 10 23 5 56 4 2 14 15 17 7 18
1 8 21 45 6 7 11 17 3 2 10 4 23 25
1 9 5 40 3 4 21 35 16 18 2 6 11 12
1 9 14 26 4 2 11 5 3 12 27 34 7 28
1 9 21 25 3 4 8 5 6 16 2 36 14 33
1 10 22 34 27 12 3 4 2 14 24 5 8 17
1 10 48 9 19 4 8 6 7 17 3 2 34 15
1 12 48 6 2 38 3 22 7 10 11 5 4 14
15 11178.60 0 なし
16 112038.24 0 なし
17 1122292.52 6 1 2 4 8 16 32 27 26 11 9 45 13 10 29 5 17 18
1 3 12 10 31 7 27 2 6 5 19 20 62 14 9 28 17
1 7 3 15 33 5 24 68 2 14 6 17 4 9 19 12 34
1 7 12 44 25 41 9 17 4 6 22 33 13 2 3 11 23
1 7 31 2 11 3 9 36 17 4 22 6 18 72 5 10 19
1 21 11 50 39 13 6 4 14 16 25 26 3 2 7 8 27

8. 計算時間

ボールの個数 $n$ と計算時間の関係を以下に示す。計測可能な範囲においては,ボールが一個増える毎に計算時間が十倍になるという綺麗な関係が成立しているように見える。

計算には ThinkPad X201(i5-540M,2.53GHz)を用いた。おそらく最新のプロセッサならこの3倍近く速いだろうが,$n = 17$ 個の場合は二週間ほどかかるため,それだけの長期間メイン PC を占有する訳にはいかず,古いサブノート PC を使わざるを得なかったのだ。

Visual C++ には 32bit コードを生成するコンパイラと 64bit 用の二種類あるが,上記の計算時間は 32bit 用を使用したものである。64bit コードのほうが使えるレジスタが倍増するため速くなるのではと思われたが,期待に反して 32bit 用のほうが早かったのだ。

C# でも組んでみたが,残念ながら C 言語ネイティブの1.5倍ほど時間がかかった。まだまだである。

9. 次回

次回へ続きます。

2
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
2
1