これは Theorem Prover Advent Calendar 2016 の5日目の記事です。
はじめに
オフラインリアルタイムどう書くというプログラマ向けイベントがあります。オフラインで集まって、主催者が出したプログラミング課題を参加者が制限時間内に解いて、回答を発表しあうというものです。そこで過去に出題された課題に対してCoq/SSReflectを使った回答を試みた、というのがこの記事の趣旨です。さすがにイベントで回答するには時間が足りなさすぎるので、自宅で解きました。この記事で扱う問題はこれです。
五角形で敷き詰められた升目を移動するためのコマンド列が与えられるので、閉路を見つけなさいという問題です。詳細はリンク先をご覧ください。回答を試みたところ、結果的にこうなりました。
- 自分なりに回答した。問題を自分なりに定式化して、その仕様を満たすことを証明した。
- しかし元の問題と微妙に異なる問題を解いていたことがわかった。修正するモチベが尽きたのでそのまま…。
- 定式化に迷った。この証明をもって健全性の保証とできるか怪しい。
なんとも塩っぱい結果ですが、本人は教科書等の手がかりなしにCoqコードを数百行程度書くのは初めてで、内容はともかく最後まで書き切り、いくらかの学びを得られたので満足しています。
環境
- Coq 8.5pl2
- SSReflect 1.6
- Ubuntu 16.04 LTS on MacOSX (Vagrant)
回答
これが自分が出した回答です。以後、おそらく上級者の皆さんにとっては冗長なコードが続きます。
実装
僕はここで示されている升目を平面に六角形を敷き詰めたものと解釈しました。1つの六角形は4つの五角形を持ち、六角形の位置をxy座標(整数)で、六角形内の五角形の位置をC0〜C3までの値で、向きをD0〜D4までの値で表しています。
Inductive cell := C0 | C1 | C2 | C3.
Inductive dir := D0 | D1 | D2 | D3 | D4.
Definition node := (cell * Z * Z)%type.
Definition state := (node * dir)%type.
状態遷移(ここでは五角形間の移動)は関数で愚直に実装しました。cell * Z * Z * dirの整数部分は元の六角形の位置との相対座標です。規則性はないのでしょうか？
Definition link := (cell * Z * Z * dir)%type.
Definition penta_graph c d : link :=
(match c with
| C0 => match d with
| D0 => (C2, 0, 0, D2) | D1 => (C3, 0, 0, D3) | D2 => (C1, -1, -1, D4)
| D3 => (C3, 0, -1, D0) | D4 => (C1, 0, 0, D1)
end
| C1 => match d with
| D0 => (C2, 0, 0, D3) | D1 => (C0, 0, 0, D4) | D2 => (C2, 0, -1, D0)
| D3 => (C3, 1, 0, D1) | D4 => (C0, 1, 1, D2)
end
| C2 => match d with
| D0 => (C1, 0, 1, D2) | D1 => (C3, 0, 0, D4) | D2 => (C0, 0, 0, D0)
| D3 => (C1, 0, 0, D0) | D4 => (C3, 1, 1, D2)
end
| C3 => match d with
| D0 => (C0, 0, 1, D3) | D1 => (C1, -1, 0, D3) | D2 => (C2, -1, -1, D4)
| D3 => (C0, 0, 0, D1) | D4 => (C2, 0, 0, D1)
end
end)%Z.
自分で定義したInductive Typeをreflectionするための補題が続きます。
Definition eqb_cell c1 c2 :=
match (c1, c2) with
| (C0, C0) => true | (C1, C1) => true | (C2, C2) => true | (C3, C3) => true | _ => false
end.
Lemma eqb_cellP : Equality.axiom eqb_cell.
Proof. rewrite /Equality.axiom. by case; case; constructor. Qed.
Canonical cell_eqMixin := EqMixin eqb_cellP.
Canonical cell_eqType := Eval hnf in EqType cell cell_eqMixin.
直進、時計回り、反時計回りの回転といったコマンドにしたがって状態遷移させる部分です。
Definition turn_dir_ccw d : dir :=
match d with
| D0 => D1 | D1 => D2 | D2 => D3 | D3 => D4 | D4 => D0
end.
Definition turn_dir_cw d : dir :=
match d with
| D0 => D4 | D1 => D0 | D2 => D1 | D3 => D2 | D4 => D3
end.
Definition move_ccw (s : state) : state :=
let: (n, d1) := s in let: d2 := turn_dir_ccw d1 in (n, d2).
Definition move_cw (s : state) : state :=
let: (n, d1) := s in let: d2 := turn_dir_cw d1 in (n, d2).
Definition move_forward (s : state) : state :=
let: ((c1, x, y), d1) := s in
let: (c2, dx, dy, d2) := penta_graph c1 d1 in
((c2, x + dx, y + dy), d2)%Z.
この問題では閉路を升目の移動距離に相当するインデックスで表現する必要があるので、コマンド列を分けて、1部分列が隣の五角形への移動に対応するようにした方が扱いやすいと考えました。1部分列を切り出す関数を split_action で定義しています。
ここでは直進、時計回り、反時計回りのコマンドをそれぞれ Fw, Tn Cw, Tn ccw と表現しました。
Inductive turn := Cw | Ccw.
Inductive command := Fw | Tn of turn.
Definition action := seq turn.
Function split_action (cs : seq command) : option (action * seq command) :=
match cs with
| [::] => None
| Fw :: cs1 => Some ([::], cs1)
| Tn t :: cs1 =>
match split_action cs1 with
| None => None
| Some (ts, cs2) => Some (t :: ts, cs2)
end
end.
これを使ってコマンド列から部分列(アクションとします)の列を作る関数を定義します。例えば、[:: Fw, Tn Cw, Fw, Tn Ccw, Tn Ccw, Fw ]というコマンド列を [:: [::], [:: Cw], [:: Ccw, Ccw] ] のように変換します。
ここで本題から離れますが、FunctionというVernacular Commandを使うと、整礎関係を使った関数定義を比較的苦しまずに実装できます。今回の実装で気づきました。整礎関係を使って関数定義するときは 標準ライブラリのRecdef モジュールをあらかじめ Import する必要があります。詳細はマニュアルを参照してください。
Lemma split_action_size cs1 cs2 a :
split_action cs1 = Some (a, cs2) -> size cs1 = size a + size cs2 + 1.
Proof.
elim: cs1 cs2 a => // c cs1 IHcs1 cs2 a H0. elim: c H0 => /=.
- case => <- -> /=. by rewrite -addn1.
- move H0: (split_action cs1) => x. case: x H0 => //.
move=> x H0 t. case: x H0 => a1 cs H0. case=> H1 H2.
move/IHcs1 :H0. rewrite -H1 H2 => /= ->. ring.
Qed.
Function convert_to_actions (cs : seq command) {measure size cs} : seq action :=
match split_action cs with
| None => [::]
| Some (a, cs') => a :: convert_to_actions cs'
end.
Proof.
move=> cs1 p a cs2 H0. move/split_action_size => ->. apply/ltP.
by rewrite [size a + size cs2]addnC -addnA -{1}[size cs2]addn0 ltn_add2l addn1 ltn0Sn.
Qed.
上記関数を使って生成されたアクションを用いた状態遷移関数です。
Function perform_action s a :=
match a with
| Cw :: a' => perform_action (move_cw s) a'
| Ccw :: a' => perform_action (move_ccw s) a'
| [::] => move_forward s
end.
後で同じ状態のインデックスを得るために状態列を作っておきます。
Fixpoint process' s aa :=
s :: match aa with
| [::] => [::]
| a :: aa' => process' (perform_action s a) aa'
end.
状態列の中に同じ状態を探索します。
Function find_same_node' ss i :=
match ss with
| [::] => None
| s :: ss' =>
match findo (compare_by_node s) ss' with
| Some j => Some (i + j + 1, i)
| None => find_same_node' ss' i.+1
end
end.
Definition find_same_node ss := find_same_node' ss 0.
ここまでの関数定義を合成して終わりです。それにしても長かった。
Definition solve' cs s n := find_same_node' (process' s (convert_to_actions cs)) n.
Definition solve cs := solve' cs initial_state 0.
証明
さて、何を証明すれば良いでしょうか。これで良いか確信できないのですが、この課題の仕様を論理プログラミング的に Inductive Type を記述することで定義してやって、元のプログラムがこの仕様を満たすことを証明すれば良いのではないかと考えました。手始めに、動く仕様を得るべく、実際に Prolog で論理プログラミングしてみました。
僕は普段Prologコードを書かないので、コードにおかしな箇所があるかもしれません。ところで、これはCoqのコードに関しても言えるのですが、字句解析の部分は省略しています。今回はアルゴリズムの証明に注力していて実装し忘れていることに直前になって気づきました…。
たまにはPrologコードを書くのも楽しいですよね。動く仕様を手に入れたところで、この中で is_answer_2, is_answer_1 に相当する部分をCoq風に書き直してみます。
Inductive is_answer_2 : state -> nat -> state -> seq command -> nat -> Prop :=
| IsAnswer21 : forall cs s1 s2 n, same_node s1 s2 -> is_answer_2 s1 n s2 cs n
| IsAnswer22 : forall cs2 cs3 s1 s2 s3 a n2 n3 n4,
next_command cs2 a cs3 -> R_perform_action s2 a s3 -> n3 = n2.+1 ->
is_answer_2 s1 n3 s3 cs3 n4 -> is_answer_2 s1 n2 s2 cs2 n4.
Inductive is_answer_1 : state -> seq command -> nat -> nat -> nat -> Prop :=
| IsAnswer11 : forall cs1 cs2 s1 s2 a n1 n2 n3, next_command cs1 a cs2 ->
R_perform_action s1 a s2 -> n2 = n1.+1 -> is_answer_2 s1 n2 s2 cs2 n3 ->
is_answer_1 s1 cs1 n1 n1 n3
| IsAnswer12 : forall cs1 cs2 s1 s2 a n1 n2 n3 n4, next_command cs1 a cs2 ->
R_perform_action s1 a s2 -> n2 = n1.+1 -> is_answer_1 s2 cs2 n2 n3 n4 ->
is_answer_1 s1 cs1 n1 n3 n4.
Inductive is_answer : seq command -> nat -> nat -> Prop :=
| IsAnswer cs m n : is_answer_1 initial_state cs 0 m n -> is_answer cs m n.
関数 solve がこの「仕様」を満たすことを証明しました。証明はコードを参照してください。
Theorem solve_spec cs m n : solve cs = Some (n, m) -> is_answer cs m n.
逆向きの証明、つまり完全性については証明していません。このアルゴリズムは閉路が複数あった場合、開始インデックスの小さい方を優先的に出力しているからです。ところが、いざこのコードをOCamlに変換して実行してみたところ、課題ページに掲載されているいくつかのテストケースが通りませんでした。どうやら出題者の意図としては、終了インデックスの小さい方を優先的に出力して欲しかったようです。このことにAdvent Calendar締め切りの直前になって気づきましたが、モチベが尽きており、修正せずにそのままブログに投稿してしまいました＞＜
グラフ自体が整合性が取れているものかどうか確認すべく、「行って戻ってきたら元の位置に戻る」を証明しました。
Theorem penta_graph_consistent s : solve' [:: Fw; Fw] s 0 = Some (2, 0).
考察のようなもの
プログラムを書いたらまず単体テストを実行すべき？
証明の前に単体テストを書くというのは、何かこう負けたような気分になります。結局、信頼性の担保は伝統的なテストに依存していることに違いないからです。しかし、僕の少ない経験から言えることは、証明に比べると単体テストの実行とデバッグははるかにコストが低いということです。証明に長い時間をかけた結果、そもそもプログラムに単純なバグがあり証明不可能であるといった状況に比べれば、単体テストを書いて実行した方が、全体的に生産性が高いのではないか。今回、比較的長いコードを書いた後で、そう振り返りました。
問題の定義にも同様のことが言えて、定式化したら、まずExampleコマンドで簡単な問題と出力からなる事実を証明できるか試すべきなのでしょう。このように。
Example is_answer_example_1 : is_answer [:: Fw; Fw] 0 2.
最初からPrologコードを書くべき？
Prologコードを書いて実行できるなら、それを最終的な成果物として採用すれば良いのではないかということです。今回に関して言えばPrologが向いていると言えそうです。しかし、一般的にPrologコードは、たとえコードそのものにバグがなくても常に現実的な時間内に終わるとは限りません。問題の仕様をPrologで記述できるが、コードが遅いといった場合に今回の方針は役に立つと言えなくもないのではないでしょうか。
話は変わりますが、Prologを持ち出してまで、いわゆる関数型プログラムの動作を保証するのが筋が良いのかという疑問もあります。僕には馴染みがないのですが、今回の問題の本質は状態遷移であることから、Alloy等の形式手法の方が今回の問題に向いているのかも知れません。勉強不足のため、動作保証に何が最適か結論を得られずにいます。
Vernacular Command "Function" は有用
再帰関数を Function を使って定義しておくと、項の書き換えに便利な場合があります。説明のためにStack Overflowから例を引用し、SSReflect用に変更したものを用います。以下のような関数 f を定義して、証明の途中で定義に置き換えたいとします。
Definition g (x:nat) := x.
Fixpoint f (x : nat) := if g x is y.+1 then (if x is z.+1 then f z else 1) else 0.
Lemma test : forall (x : nat), g x = O -> f x = O.
Proof.
move=> x.
1 subgoal, subgoal 1 (ID 7)
x : nat
============================
g x = 0 -> f x = 0
ここで rewrite /f すると、 fix が現れてその後の証明が難しくなります。
1 subgoal, subgoal 1 (ID 8)
x : nat
============================
g x = 0 ->
(fix f (x0 : nat) : nat :=
match g x0 with
| 0 => 0
| _.+1 => match x0 with
| 0 => 1
| z.+1 => f z
end
end) x = 0
そこで関数 f を Function を使って定義します。
Function f (x : nat) := if g x is y.+1 then (if x is z.+1 then f z else 1) else 0.
すると、f の定義と共に6つの補題が自動生成されます。
f is defined
f is recursively defined (decreasing on 1st argument)
f_equation is defined
f_ind is defined
f_rec is defined
f_rect is defined
R_f_correct is defined
R_f_complete is defined
このうち f_equation は以下のような補題です。
f_equation
: forall x : nat,
f x = match g x with
| 0 => 0
| _.+1 => match x with
| 0 => 1
| z.+1 => f z
end
end
rewrite f_equation すると、素直に定義に置き換えてくれます。
1 subgoal, subgoal 1 (ID 255)
x : nat
============================
g x = 0 -> match g x with
| 0 => 0
| _.+1 => match x with
| 0 => 1
| z.+1 => f z
end
end = 0
後は rewrite /g => -> //. で無事終了。
Function は書き換えの他にも関数定義から健全性と完全性の定理を自動生成してくれることがあります。こちらはどうやら上手くいく場合といかない場合があるようです。
まとめ
- プログラム課題をCoqで実装し、問題を定式化し、その仕様を満たすことを証明した。
- しかしそれをもって健全性の保証として良いか疑問が残る。
- 著者はプログラムを書いたら単体テストを実行すべきではと考えている。
- Prologは有用。でも常に実行できるとは限らないから、その時のために、いわゆる関数型プログラムとの同値性を保証できた方が良いと著者は考えている。
- Vernacular Command
Functionは有用。積極的に利用しよう。
ここまで読んでいただき本当にありがとうございます。