ナンプレを線形計画法で解く際のロジックについて

Question

こちらのサイトを参考にナンバープレイス（通称ナンプレ）を線形計画法で解こうと考えています。
まず整数計画法のモデルを作るためにナンプレのルールを以下とします。

ルール
1. 横一列に同じ数字はこない
2. 縦一列に同じ数字はこない
3. 3x3 のブロック内に同じ数字はこない

モデルに関しては

0-1 整数計画
目的関数は不要
としています

import pulp
from itertools import product
from pprint import pprint

def main():
    # 問題表(0が空白マス)
    Board = [[0, 0, 0, 0, 4, 0, 6, 0, 0],
             [7, 5, 0, 0, 0, 0, 0, 0, 0],
             [0, 0, 0, 0, 0, 0, 0, 0, 0],
             [0, 0, 0, 0, 0, 5, 0, 0, 0],
             [0, 0, 9, 0, 0, 1, 0, 3, 0],
             [0, 0, 6, 0, 0, 0, 4, 0, 0],
             [0, 0, 0, 0, 0, 0, 0, 1, 0],
             [0, 0, 0, 0, 3, 7, 0, 0, 5],
             [2, 0, 4, 0, 6, 0, 0, 0, 0]]


    # 問題の宣言
    problem = pulp.LpProblem()

    # 定数の宣言
    line = list(range(9))
    
    # 変数の宣言
    x = [[[pulp.LpVariable(f'x{i}{j}{k}', cat='Binary')
           for k in line] for j in line] for i in line]
    
    # 目的関数
    # 今回はなし

    # 制約条件
    # すでに埋まっているマスの情報を入力
    for i, j in product(line, line):
        
        if Board[i][j] > 0:
            
            problem += x[i][j][Board[i][j]-1] == 1
    
    # 数独のルールに合うように制約を追加
    for j, k in product(line, line):
        #print(i,j,k)
        problem += pulp.lpSum(x[i][j][k] for i in line) == 1

    for k, i in product(line, line):
        print(i,j,k)
        problem += pulp.lpSum(x[i][j][k] for j in line) == 1

    for i, j in product(line, line):
        #print(i,j,k)
        problem += pulp.lpSum(x[i][j][k] for k in line) == 1

    three_block = [[ii, ii+1, ii+2] for ii in [0, 3, 6]]
    for a, b in product(three_block, three_block):
        for k in line:
            problem += pulp.lpSum(x[i][j][k] for i in a for j in b) == 1

    # ソルバーの実行
    # problem.solve(pulp.GUROBI())
    problem.solve()
    
    
    # 結果の表示
    print('Result')
    for i, j, k in product(line, line, line):
        if x[i][j][k].value() == 1:
            
            Board[i][j] = k + 1

    pprint(Board)

if __name__ == '__main__':
    main()

参考

わからないこと

ライブラリに投げているのでproblem.solve()の中でどのような動きが起きているかわからないです

自分の考えとしては
貪欲法を用いて総当たりで最悪9!*9通りを総当たりして枝狩り？を行う。途中で行けそうでないなら
やめればもう少し少ない計算量で行けるのかなと思いました。(そもそも総当たりが9!*9通りであってるのかな・・・)
どうやって最適解を見つけているのでしょうか

Comments

[tiitoi]

maisumakun さんのご回答の通り、実用上は自前実装ではなく、ソルバーを使うのが推奨されます。
ソルバーは高度に最適化されているのですが、元となるアルゴリズムが気になるのであれば、「分枝限定法」で調べてみるとよいでしょう。

[eggpol]

なるほど、「分技限定法」を調べてみます。ありがとうございました。

Answer 1

ライブラリに投げているのでproblem.solve()の中でどのような動きが起きているかわからないです

それを考えなくて済むのがソルバーに任せるメリットだと思うのですが、どうでしょうか？

Answer 2

調べるさいの道標として回答します．（CBC のソースコードを見ていないので間違っていたらすみません）．

混合整数計画ソルバーとして PuLP から Cbc (Coin-or branch and cut) を呼び出しているので problem.solve() でどのようなアルゴリズムが使用されているのかは Cbc を調べる必要があります．CBC User Guide には 分枝切除法（branch and cut method） をメインで使用していると書いているので「どのような動きをしているのか」を知る手がかりとして 分枝切除法 を調べるのをオススメします．ただし，一般に混合整数計画ソルバーにはいくつかのヒューリスティック（精度保証の無いアルゴリズム）を合わせて使用している場合があるので，problem.solve() を詳細に調べるときにはソースコードと CBC User Guide 等を参照する必要があると思います．