なぜ数独なのか
前回紹介した Python の線形計画法ライブラリ PuLP
http://qiita.com/mzmttks/items/82ea3a51e4dbea8fbc17
で、数独が解けるので解いてみる。
PyConJP のトーク 数理最適化によるパズルの解法
https://pycon.jp/2014/schedule/presentation/23/
や、ドキュメント
https://pythonhosted.org/PuLP/CaseStudies/a_sudoku_problem.html
に書かれているので自分でもやってみた。
という話。
数独って何?
1-9 までの数字を制約に基づいて埋めるパズルゲーム。
これらが制約
- 横一列に同じ数字はこない
- 縦一列に同じ数字はこない
- 3x3 のブロック内に同じ数字はこない
ドキュメントの例
数独
これを解くと、
数独解
こうなる。で、これが出れば OK
数独を線形計画問題に変換する
いわゆるモデル化。これがこの手の問題のキモで、結構感動した。
まずは PuLP で問題オブジェクトを作る
import pulp
prob = pulp.LpProblem('Sudoku', pulp.LpMinimize)
prob += 0 # No objective function
目的関数は特にないので、定数の 0 とか設定しておく。
変数の設計
キーアイデアは、 0 か 1 の値のみとるバイナリ変数を大量に使う ということ。
具体的には、横9マス, 縦9マス, 入る値が 1-9 の 9 種類なので、
9 x 9 x 9 = 729 個の変数を用意して、
横 x, 縦 y の位置のマスの数字が v のとき変数 Cell_v_x_y = 1
そうでないとき Cell_v_x_y = 0
とする。
つまり、1つのマスごとに 9個の変数を用意しておいて、
それらのうち 1個だけが 1, それ以外の 8個が 0 になる。
変数たちを生成するコードは下記のとおり
numbers = range(1, 10)
xs = range(1, 10)
ys = range(1, 10)
choices = pulp.LpVariable.dicts("Cell",(xs, ys, numbers),0,1,pulp.LpInteger)
これで、Cell_1_1_1, ..., Cell_9_9_9 の 729 個の変数たちが生成されて、
choices[1][1][1]
のような感じでアクセスできるようになった。
変数は整数で、0か1の値をとる。
ゲームのルールから制約条件の設計
1 つのマスに入る値は 1 つだけ
たとえば、座標 (2, 3) のマスに入るのは 1-9 のうち 1つだけなので、
このような式になる。
Cell_1_2_3 + Cell_2_2_3 + Cell_3_2_3 +
Cell_4_2_3 + Cell_5_2_3 + Cell_6_2_3 +
Cell_7_2_3 + Cell_8_2_3 + Cell_9_2_3 = 1
コードで書くと
ys = range(1, 10)
xs = range(1, 10)
numbers = range(1, 10)
for y in ys: # 各列と
for x in xs: # 各行について、
prob += pulp.lpSum([choices[v][x][y] for v in numbers]) == 1
# すべての数字を足すと 1 になる。
ここで、lpSum
は、配列の変数をすべて足した式を表す。
総和を表すときにシンプルになるので多用される。
既に数字が埋められているときは必ずその数字のみがマスに入る。
board[y][x]
に、既に埋められているマスの数字が入っていて、
空だったら 0 が入っているとする。
for x in range(width): # 各列と
for y in range(height): # 各列を調べて、
if board[y][x] > 0: # もし数字が 0 より大きかったら、
prob += choices[board[y][x]][x+1][y+1] == 1
# その場所の数字が 1 になるような制約を追加する。
すでに、同じマスで 1 になるのは 1 つの値だけという制約を入れたので、
値が 1 のときのことだけ考えればよい。
同じ列・行・ボックスには同じ数字は入らない
行の制約
例:1列目(y)で、数字(v)が 1 になるのは、すべての行のうちどれか1つだけ。
つまり、それらを全部足すと必ず 1 になる。
Cell_v_x_y
の変数で v=1, y=1 のとき、下記のような制約になる
Cell_1_1_1 + Cell_1_2_1 + Cell_1_3_1 +
Cell_1_4_1 + Cell_1_5_1 + Cell_1_6_1 +
Cell_1_7_1 + Cell_1_8_1 + Cell_1_9_1 = 1
for v in numbers:
for y in ys:
prob += pulp.lpSum([choices[v][x][y] for x in xs]) == 1
列とボックスは省略。
コード
import pulp
import numpy
# make problem
board = """
530070000
600195000
098000060
800060003
400803001
700020006
060000280
000419005
000080079"""
# board = map(lambda line: map(int, line.rstrip()), open("sudoku.txt").readlines())
board = [map(int, list(line)) for line in board.strip().split("\n")]
print board
width = len(board[0])
height = len(board)
## initialize
# Create a problem
prob = pulp.LpProblem('Sudoku', pulp.LpMinimize) # or minimize
## objective
# No objective function
prob += 0
# make
numbers = range(1, 10)
xs = range(1, 10)
ys = range(1, 10)
choices = pulp.LpVariable.dicts("Choice",(xs, ys, numbers),0,1,pulp.LpInteger)
## constraints
# Add given-value constraints
for x in range(width):
for y in range(height):
if board[y][x] > 0:
prob += choices[board[y][x]][x+1][y+1] == 1
# one cell must have one value
for y in ys:
for x in xs:
prob += pulp.lpSum([choices[v][x][y] for v in numbers]) == 1, ""
# horizontal line must have different values
for v in numbers:
for y in ys:
prob += pulp.lpSum([choices[v][x][y] for x in xs]) == 1
for x in xs:
prob += pulp.lpSum([choices[v][x][y] for y in ys]) == 1
for x in [1, 4, 7]:
vs = []
for y in [1, 4, 7]:
print [[choices[v][x+i][y+j] for i in range(3) for j in range(3)]]
prob += pulp.lpSum([[choices[v][x+i][y+j] for i in range(3) for j in range(3)]]) == 1
# print prob
s = prob.solve()
# solve it
print pulp.LpStatus[s]
# print result
for y in ys:
for x in xs:
for v in numbers:
if choices[v][x][y].value() == 1:
print v,
print
実行結果は こんな感じになって、解けましたとさ。
5 3 4 6 7 8 9 1 2
6 7 2 1 9 5 3 4 8
1 9 8 3 4 2 5 6 7
8 5 9 7 6 1 4 2 3
4 2 6 8 5 3 7 9 1
7 1 3 9 2 4 8 5 6
9 6 1 5 3 7 2 8 4
2 8 7 4 1 9 6 3 5
3 4 5 2 8 6 1 7 9