数独の平凡な解法(Mathematica)

数独の平凡な解法(C言語)という記事を書いたら、「わかりやすさを重視する時はC言語なのですか」と言われてしまったので、いつものようにMathematicaでも書いておきます。

canBePlaced[board_, {row_, col_}, b_, x_] := And[
  Not[MemberQ[board[[row]], x]]  ,
  Not[MemberQ[board[[All, col]], x]] ,
  With[{
    r = 1 + b Quotient[row, b, 1],
    c = 1 + b Quotient[col, b, 1]},
   Not[MemberQ[board[[r ;; r + b - 1, c ;; c + b - 1]], x, 2]]]]

check[board_] := With[{next = Position[board, 0]}, 
  If[Length@next == 0, Sow@TableForm@board, (* solution is found *)
   Map[If[canBePlaced[board, First@next, Sqrt@Length@board, #], 
      check[ReplacePart[board, First@next -> #]]] &, Range@Length@board]]]

実行例は以下の通り。

Reap[check@{
   {1, 0, 0, 0, 0, 7, 0, 9, 0},
   {0, 3, 0, 0, 2, 0, 0, 0, 8},
   {0, 0, 9, 6, 0, 0, 5, 0, 0},
   {0, 0, 5, 3, 0, 0, 9, 0, 0},
   {0, 1, 0, 0, 8, 0, 0, 0, 2},
   {6, 0, 0, 0, 0, 4, 0, 0, 0},
   {3, 0, 0, 0, 0, 0, 0, 1, 0},
   {0, 4, 0, 0, 0, 0, 0, 0, 7},
   {0, 0, 7, 0, 0, 0, 3, 0, 0}}][[2]]

  1 6 2 8 5 7 4 9 3
  5 3 4 1 2 9 6 7 8
  7 8 9 6 4 3 5 2 1
  4 7 5 3 1 2 9 8 6
{{9 1 3 5 8 6 7 4 2}}
  6 2 8 7 9 4 1 3 5
  3 5 6 4 7 8 2 1 9
  2 4 1 9 3 5 8 6 7
  8 9 7 2 6 1 3 5 4

数独で見るRuby(とMathematica)のパワーと表現力で紹介した方法のほうが、他の問題に転用しやすい汎用的なものだと思いますが、数独が解ければそれで十分という場合には、この方法のほうがわかりやすくていいでしょう。

数独の平凡な解法(C言語)

数独の解き方をいくつか紹介したことがあります。

困ったことに、「どれも変態でわかりにくい、ふざけんな」という批判を受けたので、みんな大好きC言語で、平凡に書いてみました(実行結果)。

数独の平凡な解法(Mathematica)

世界一難しい数独

解けたら天才! フィンランドの数学者が作った「世界一難しい数独」というのがあるらしい。

こちらの問題、実はフィンランドの数学者Arto Inkala氏が作成した問題で、曰く「世界でもっとも難しい数独」とのこと

1847534511この問題の作者であるArto Inkalaさんは、これまで「史上最難」と言われいた問題「AI Escargot」の作者でもある。その彼が言うのだから、間違いないと言いたいところだが、そもそも、数独の難しさに客観的な指標なんてないはず。合理的なものに限っても、たまたま知っていた定石で一発、なんてこともあり得る。だから、「ある解き方においては最も難しい」という限定的なことしか言えないのだ。彼自身、著書“AI Escargot”の中で次のように述べている。

it is often impossible to clearly identify the clearly most difficult puzzle in the world.

それにも拘わらず「世界でもっとも難しい」なんて言うことの大人の事情は放っておいて、せっかくだから手元にあるアルゴリズム(単純なしらみつぶし)でどうなるかを見てみよう。解を見つけるまでにどれだけの候補をチェックしなければならないかを調べる。Mathematica版なら全体をTraceで囲んで、結果の長さを見ればよい。

problem={
  {0, 0, 5, 3, 0, 0, 0, 0, 0},
  {8, 0, 0, 0, 0, 0, 0, 2, 0},
  {0, 7, 0, 0, 1, 0, 5, 0, 0},
  {4, 0, 0, 0, 0, 5, 3, 0, 0},
  {0, 1, 0, 0, 7, 0, 0, 0, 6},
  {0, 0, 3, 2, 0, 0, 0, 8, 0},
  {0, 6, 0, 5, 0, 0, 0, 0, 9},
  {0, 0, 4, 0, 0, 0, 0, 3, 0},
  {0, 0, 0, 0, 0, 9, 7, 0, 0}};

Reap@Length[Trace[search[{problem}, Join[#2, #1] &, False]]]

{80070, {{1 4 5 3 2 7 6 9 8}}}
          8 3 9 6 5 4 1 2 7
          6 7 2 9 1 8 5 4 3
          4 9 6 1 8 5 3 7 2
          2 1 8 4 7 3 9 5 6
          7 5 3 2 9 6 4 8 1
          3 6 7 5 4 2 8 1 9
          9 8 4 7 6 1 2 3 5
          5 2 1 8 3 9 7 6 4

となるので、この問題をMathematicaで解くときの、Traceの長さは80070。

一方、AI Escargotはと言うと、

problem={
  {1, 0, 0, 0, 0, 7, 0, 9, 0},
  {0, 3, 0, 0, 2, 0, 0, 0, 8},
  {0, 0, 9, 6, 0, 0, 5, 0, 0},
  {0, 0, 5, 3, 0, 0, 9, 0, 0},
  {0, 1, 0, 0, 8, 0, 0, 0, 2},
  {6, 0, 0, 0, 0, 4, 0, 0, 0},
  {3, 0, 0, 0, 0, 0, 0, 1, 0},
  {0, 4, 0, 0, 0, 0, 0, 0, 7},
  {0, 0, 7, 0, 0, 0, 3, 0, 0}};

Reap@Length[Trace[search[{problem}, Join[#2, #1] &, False]]]

{71766, {{1 6 2 8 5 7 4 9 3}}}
          5 3 4 1 2 9 6 7 8
          7 8 9 6 4 3 5 2 1
          4 7 5 3 1 2 9 8 6
          9 1 3 5 8 6 7 4 2
          6 2 8 7 9 4 1 3 5
          3 5 6 4 7 8 2 1 9
          2 4 1 9 3 5 8 6 7
          8 9 7 2 6 1 3 5 4

となる。解くときのTraceの長さ71766は、先の値80070より小さい。ということは、ここで採用しているアルゴリズムを基準にすれば、今回発表された問題は、確かにAI Escargotより難しい。

しかし、このアルゴリズムを基準にするなら、もっと難しい問題はある(Mathematicaでは時間がかかるからC++版で試す。オプション「-pg」を付けてコンパイルして実行し、gprofでプロファイルを見ればよい)。MathematicaでもC++でも、わざわざコードを書き換えたりしないのがポイント。

要するに、どの隅から解いていくかという話なのだが、ここで使っている方法の場合、どちらの問題も、右下から左下に向かって解くときに最も難しくなる。再帰的な関数searchが、今回発表された問題180度回転したものの場合は110874回、AI Escargot180度回転したものの場合は1502267回呼ばれる(どちらもCodepadでは解けない)。

つまり、ここで採用した基準では、世界一難しいのは今回発表された問題ではなく、依然としてAI Escargot(を180度回転したもの)だということになる。

追記:どちらも自力で解けた友人によると、AI Escargotの方が難しかったそうです。もう一つの方は30分ほどだったとか。

SQLによる数独の解法 #deim2010

2/28から3/2にかけて開催されたDEIM2010で発表してきました。

ワークショップの様子

  1. Togetter – まとめ「DEIM2010 初日」
  2. Togetter – まとめ「DEIM2010 2日目」
  3. Togetter – まとめ「DEIM2010 最終日」

「データベースの授業で使いたい」というコメントを多くの方からいただいたので、そのためのマテリアルを公開します。

論文の最終版は後で公開します。

論文:矢吹太朗, 佐久田博司. SQLによる数独の解法とクエリオプティマイザの有効性. 日本データベース学会論文誌. Vol.9, No.2, pp.13–18, 2010.

参考:数独の平凡な解法(C言語Mathematica

プログラマの道具箱(深さ優先探索と幅優先探索) C++編

参考:数独の平凡な解法(C言語Mathematica

Mathematicaで実装した深さ優先探索と幅優先探索をC++に移植してみましょう。

深さ優先探索と幅優先探索

探索木のノードを表現する型をSTATUS、そのリストをFRINGE、FRINGEにノードを追加するメンバ関数をCOMBINERとします。

データ構造はツリーではなくリストを使います。候補ノードを、スタックで管理すれば深さ優先探索、キューで管理すれば幅優先探索になるわけですが、std::stackとstd::queueにおいて、要素を取り出すメソッド名が同じではないため多態性を利用できません。ですから、ノードはstd::listで管理し、要素はメンバ関数pop_front()で常に先頭から取り出すことにします。

要素を追加する際に、メンバ関数push_front()を使えば深さ優先探索に、push_back()を使えば幅優先探索になります。探索に使う関数search()は、引数としてこれらのメンバ関数へのポインタを取るようにしましょう。

道具箱には次のようなコードを入れておけばいいでしょう。

#include <iostream>
#include <list>

using namespace std;

typedef ... STATUS; //ノードのための型
typedef list<STATUS*> FRINGE; //ノードを管理するオブジェクトの型
typedef void (FRINGE::*COMBINER)(FRINGE::const_reference); //メンバ関数へのポインタ

//解を報告する手続き
void report(STATUS* s)
{
  ...
}

//解かどうかを判定する述語
bool goal(STATUS* s)
{
  ...
}

//探索を深める
void expand(FRINGE* fringe, STATUS* s, COMBINER combiner)
{
  ...
  STATUS* newStatus=new STATUS...
  ...
  (fringe->*combiner)(newStatus);
}

//探索の本体
void search(FRINGE* fringe, COMBINER combiner, bool findAll)
{
  if (fringe->size()!=0) {
    STATUS* s=fringe->front();
    fringe->pop_front();
    if (goal(s)) {
      report(s);
      if (!findAll) {
        while (!fringe->empty()) {
          STATUS* tmp=fringe->front();
          fringe->pop_front();
          delete tmp;
        }
      }
    }
    else expand(fringe, s, combiner);
    delete s;
    search(fringe, combiner, findAll);
  }
}

int main ()
{
  STATUS* s=new STATUS...
  ...
  FRINGE* fringeP=new FRINGE();
  fringeP->push_back(s);
  search(fringeP, &FRINGE::push_front, true); //depth-first search
  //search(fringeP, &FRINGE::push_back, true); //breadth-first search
  return 0;
}

数独を解く場合には、typedef int STATUS;でいいでしょう。関数については、実行結果を参照してください(実行する場合には、スタックサイズを大きくする必要があるかもしれません)。

コメントと空行を取り除くと、(Mathematicaの時は含めなかった)実行部分を含めて約100行です。

C++に関する必要な知識の大部分は、拙著『Microsoft Visual C++入門』にまとめてありますが、深さ優先・幅優先を指定するために必要なメンバ関数へのポインタの使い方は確認しておきましょう。

クラスFooのオブジェクトのメンバ関数へのポインタには、次のように名前(ここではBAR)を付けておくのが簡単です。

typedef 戻り値の型 (Foo::*BAR)(パラメータリスト)

Fooののメンバ関数funcへのポインタbarは、次のように取得できます。

BAR bar=&Foo::func;

次のようにして、クラスFooのオブジェクトfooのメンバ関数を、ポインタbarを使って呼び出します。

(foo->*bar)(引数リスト)

参考:プログラマの道具箱(深さ優先探索と幅優先探索) Mathematica編