AtCoder Beginner Contest の E 問題で2回続けてダイクストラ法が出題される昨今,いかがお過ごしでしょうか?1
その2つの E 問題のうち,先に出題された E - Come Back Quickly はいかにも典型的だなあといった問題なので良いのですが,後に出題された E - Train は一筋縄ではいかないものでした.典型的な最短経路問題は辺の長さが定数であるようなグラフを対象としていますが,その問題では辺の長さが,その辺に辿り着くまでの経路長に依存するようなグラフに対してダイクストラ法を適用するのが想定解だったからです.
本記事では,頂点に辿り着くまでの経路長 に頂点 から頂点 に伸びる辺の長さを加えたものが関数 として与えられた際に, が単調増加関数かつ を満たせばダイクストラ法が最短経路を返すことを解説します.加えて,ナイーブな実装ではありますが,そのようなグラフに対してのダイクストラ法が最短距離を返すことを Coq で証明します.
ダイクストラ法のお気持ち
辺の長さが変化するようなグラフへダイクストラ法を拡張する前に,まずは普通のグラフを対象としたダイクストラ法がなぜ最短経路を求められるか振り返ってみましょう.最短経路問題を解くアルゴリズムたらしめる性質を詳かにすれば,それを一般化して辺の長さが変化するようなグラフにも適用できることでしょう.
一般に知られているダイクストラ法とは,長さが負の辺を持たない有向グラフと,始点となる頂点が与えられた時に,全ての頂点について始点からの最短距離2を求めるアルゴリズムです.大まかな実装としては以下のようになります.
- 始点への距離を0,それ以外の頂点への距離を無限大で初期化する
- 最短距離がまだ確定していない頂点の集合 を全体集合で初期化する
- が空集合であれば,全ての頂点への最短距離が確定したので終了する
- に属する頂点のうち,距離が最小となる頂点 は最短距離が確定したので, から取り除く
- に隣接する全ての頂点 について,元々記録していた への距離と, を経由して に向かう際の距離を比較し,小さい方を への距離とする
- 3 に戻る
参考のために,OCaml によるナイーブな3実装を示します.
let dijkstra (* 頂点数 n *) (n : int) (* 頂点 u ∈ [0, n) と頂点 v ∈ [0, n) を受け取って,u から v へ向かう辺の長さを返す関数 *) (adj : int -> int -> int) (* 始点 s ∈ [0, n) *) (s : int) (* 終点 t ∈ [0, n) を受け取って,s から t への最短距離を返す関数 *) : (int -> int) = (* ダイクストラ法のループ *) let rec dijkstra_rec (* 最短距離が確定していない頂点のリスト *) (vs : int list) (* 頂点 v ∈ [0, n) を受け取って,現在の v への距離を返す関数 *) (d : int -> int) (* 終点 t ∈ [0, n) を受け取って,s から t への最短距離を返す関数 *) : (int -> int) = match vs with (* 全ての頂点への最短距離が確定したので,今覚えているのが最短距離 *) | [] -> d (* vs が少なくとも頂点 v を含む場合 *) | v :: _ -> (* リスト vs を線形探索し,距離が最小となる頂点 u を求める *) let u = List.fold_left (fun u v -> if d u <= d v then u else v) v vs in dijkstra_rec (* u を vs から取り除く *) (List.filter (( <> ) u) vs) (* 全ての頂点 v について,元々の v への距離と u を通って v に向かう距離を比較し,小さい方で更新 *) (fun v -> min (d v) (d u + adj u v)) in dijkstra_rec (* 最短距離が確定していない頂点のリストを,[0; 1; 2; ... ; n - 1] で初期化 *) (List.init n Fun.id) (* 始点への距離を 0,それ以外の距離を最大の整数で初期化 *) (fun v -> if v = s then 0 else max_int)
では,なぜ長さが負の辺を持たない有向グラフは,ダイクストラ法で最短経路問題を解くことができるのでしょうか?
まず気になるのは,何故最短距離が確定していない頂点の集合 の要素のうち,距離が最小となる頂点 は最短距離が確定したと見做して良いのかですが,これは対象となるグラフは全ての辺の長さが非負であることから正当化されます.全ての辺の長さが非負なので,一度 から他の頂点へ移動して紆余曲折あってまた に戻ってくるような経路の距離は,最初の までの距離より必ず大きくなりますし, に属する別の頂点 から紆余曲折あって に至る経路は,元々 のものより長かった の経路長よりさらに長くなるので,そのような経路は想定しなくて良いという訳です.
次に気になるのは,何故頂点 への最短距離が確定したタイミングでだけ を通って他の頂点 に向かう経路を考慮すれば十分かについてですが,これは から に向かう辺の長さ を加算する関数 が単調増加する ことから正当化されます. への距離は検討した経路長のうち最も小さいものを残していくので実行が進むごとに減少していくのですが,まだ最短だと決まっていない段階の距離 で頂点 に向かったところで,その経路長 は最短だと決まった後の距離 で に向かう経路長 より短くなることはない ので,そのような経路は考慮しなくて良い訳ですね.
これらの工夫による時間計算量の削減たるや凄まじく,負の辺も考慮しなくてはならないベルマンフォード法は を頂点数 を辺の数とすると 回も経路長の比較を行うのに対して,ダイクストラ法は 回の比較で済んでいます.この効率的なアルゴリズムを色んな場面で使えると嬉しいと思いませんか?
辺のコストが変化するグラフへの拡張
辺の長さが負ではない定数で与えられるグラフに対してダイクストラ法が最短経路を与える理由を確認したところで,それを一般化して,どのような性質が成り立てば辺のコストが変化してもダイクストラ法を適用できるか考えてみましょう.
まず頂点 から頂点 に向かう辺の長さを, に至るまでの経路長 についての関数 と置い……ても良いのですが,実際にプログラミング言語で実装する時のことを考えて, から に向かう辺の長さを までの経路長 に加えたものを と置きます.経路長が float
だったり色んな型の場合も考慮するとなると加算の実装もダイクストラ法の実装に与えたくなるんですが,関数2つ渡すよりは だけで全部やってくれると嬉しいので.
全ての辺の長さが非負であることを についての条件で書くと,(全ての頂点 ,,経路長 について) になります. が から に向かう辺の長さなら と書けたんですが,実際はそれに を加えたものを と定義しているので,等価な不等式にするために左辺にも を加えて にしなくてはなりません.
また,辺の長さが定数であった時は加算の性質から自然に成り立ったのですが,頂点 への最短距離が確定したタイミングでだけ を通って他の頂点 に向かう経路を考慮すれば十分なことを言うために, は単調増加関数でなくてはなりません.もし が単調増加関数ではない,つまり なのに を満たすようなことがあれば, にとって最短ではない経路の後に辺 を通ったものが にとっての最短経路だった,なんてことが起こりかねません.ダイクストラ法は一種の貪欲法なので,部分構造最適性を満たして欲しいんですよね.
まとめると,頂点 から頂点 に伸びる辺の長さが に辿り着くまでの経路長 の関数になる場合,それらの和 が満たすべき性質は以下の通りです.
- は単調増加関数
とはいえ,こうして挙げた条件で十分なのかは結構気になるものです.こういう時に僕がお勧めしたいのは,Coq で証明を書き下して,洞察が正しいことを機械に保証してもらうことです.
Coq による検証
それまでの経路長によって辺の長さが変化するグラフにダイクストラ法を適用できる条件が分かったところで,本当にその条件を満たせばダイクストラ法で最短経路を求められることを Coq で検証しましょう.
まず,対象となるグラフを定義します.頂点の全体集合 V
は有限集合で,経路長の集合 C
は全順序集合です.
From mathcomp Require Import all_ssreflect. Section Dijkstra. Import Order.TTheory. Variable V : finType. Variable C : orderType tt.
グラフに含まれる辺は隣接行列のノリで定義します.adj v u
は,それまでの経路長 c
を受け取って,頂点 u
から 頂点 v
に伸びる辺の長さを加えたものを返す関数で,前章での に相当しています.引数の並びに違和感があるかもしれないですが,adj x w (adj w v (adj v u c))
みたいな感じに続けて適用するのを考えると,始点を後にした方が違和感が少なかったのでこのような定義にしました.
Variable adj : V -> V -> C -> C.
前章で考察した,辺のコストを加算する関数 adj v u
の満たすべき性質を定義しましょう.c <= adj v u c
が成り立つことと,adj v u
が単調増加関数であることでしたよね.
Hypothesis adj_increase : forall v u c, (c <= adj v u c)%O. Hypothesis adj_monotone : forall c c', (c <= c')%O -> forall v u, (adj v u c <= adj v u c')%O.
次に,始点までの経路長を c0
として,途中通った頂点のリスト p
を受け取って終点 t
までの経路長を返す関数 cost_of_path
を定義しましょう.cost_of_path t [:: vn, ..., v2, v1, s]
は経路 の長さを表しています.経路が存在することを帰納的に定義しても良いのですが,ダイクストラ法の正当性を証明する際に最短距離が確定した頂点だけを通る経路みたいなことを表現したくなったので,こんな定義にしています.
Section WithC0. Variable c0 : C. Fixpoint cost_of_path t p := if p is v :: p then adj t v (cost_of_path v p) else c0.
この経路長の定義を用いると,与えられた経路 が頂点 から 頂点 への最短経路であることを表す述語 shortest_path t [:: vn, ..., v2, v1, s]
は次のように表現できます.last x xs
はリスト x :: xs
の最後の要素を返す関数で,last t p' = last t p
の部分は経路 p'
の始点が経路 p
と同じであることを仮定していることに注意してください.
Definition shortest_path t p := forall p', last t p' = last t p -> (cost_of_path t p <= cost_of_path t p')%O.
経路長を求める関数 cost_of_path
を用いると,ダイクストラ法によって始点 から任意の頂点への最短経路を全て求める関数 dijkstra s
は次のように定義できます.ここで,最短経路の長さではなく最短経路を求める定義にしているのは,関数が返した経路長に対応する経路が存在して最短経路になっている,みたいなことを証明するよりも,関数が返した経路が最短経路になっている,みたいなことを証明する方が楽だったからです.
Fixpoint dijkstra_rec n vs ps := if n is n.+1 then if vs is v :: _ then let u := [arg min_(u < v in vs) cost_of_path u (ps u)]%O in dijkstra_rec n (rem u vs) (fun v => if (cost_of_path v (ps v) <= adj v u (cost_of_path u (ps u)))%O then ps v else u :: ps u) else ps else ps. Definition dijkstra s := dijkstra_rec (#|V|.-1) (rem s (enum V)) (fun v => if v == s then [::] else [:: s]).
ここでダイクストラ法の初期化後のループ dijkstra_rec n vs ps
の各引数が何を表しているのか補足しておくと,n
は最短経路が確定していない頂点の数,vs
は最短経路が確定していない頂点のリスト,ps
は頂点 を受け取ってこれまでに見つかった から への最も短い経路を返す関数です.vs
があったら n
要らなくないか?って気持ちになるんですが,Coq の停止性判定はあまり賢くないので,こういう余計な引数が必要になるのはありがちです.
あと SSReflect 4 とかいうライブラリの機能を使っているので補足しておくと,rem u vs
はリスト vs
の中で最初に現れた u
を取り除いたリストを返す関数で,[arg min_(u < v in vs) cost_of_path u (ps u)]%O
はリストv :: vs
の中で cost_of_path u (ps u)
が最小となる要素 u
を返す関数,#|V|
は有限集合 V
の要素数,enum V
は有限集合 V
に含まれる要素のリストです.
ダイクストラ法の実装が定義できたので,実際にその正しさを証明していきましょう.
我々が証明したいのはダイクストラ法の実装 dijkstra s
が始点 からの最短経路を返すことなので,主定理のステートメントは次のようになります.
Theorem dijkstra_correct s : (forall v, last v (dijkstra s v) = s) /\ (forall v, shortest_path v (dijkstra c0 s v)).
もっとも,このまま帰納法による証明を試みても帰納法の仮定が弱すぎて証明が回らないので,一般の n
,vs
,ps
について dijkstra_rec n vs ps
の満たすべき性質を証明し,その系として主定理を証明します.帰納法あるあるですね.
主定理 dijkstra_correct
を含み,帰納法が回る程度に一般的な,ダイクストラ法のループ dijkstra_rec
についての補題として僕が選んだのは,次のようなステートメントです.
Lemma dijkstra_rec_correct c0 s : forall n vs, size vs = n -> uniq vs -> forall ps, (forall u, last u (ps u) = s) -> (forall v, { in ps v, forall u, u \notin vs }) -> (forall u, u \notin vs -> { in vs, forall v p, last v p = s -> (cost_of_path c0 u (ps u) <= cost_of_path c0 v p)%O }) -> (forall v p, { in p, forall u, u \notin vs } \/ v \notin vs -> last v p = s -> (cost_of_path c0 v (ps v) <= cost_of_path c0 v p)%O) -> (forall v, last v (dijkstra_rec c0 n vs ps v) = s) /\ (forall v, shortest_path c0 v (dijkstra_rec c0 n vs ps v)).
まず,n
は最短経路が確定していない頂点のリスト vs
を集合として見た時の要素数なので,size vs = n
です.vs
が重複した要素を持っているとリストの長さが集合として見た時の要素数と一致しなくなるので,各要素が一意であることも仮定します.
size vs = n -> uniq vs ->
次に,頂点 v
へのこれまでに見つかった中で最も短い経路 ps v
は,必ず始点 s
から始まる経路でなければなりません.この仮定を書かなくても良いような,何か上手い定式化は無いものでしょうか.
(forall v, last v (ps v) = s) ->
端点への最短経路が確定した時しか隣接した頂点の経路を更新しないので,頂点 v
に至るまでの経路である ps v
は,当然最短経路が確定した頂点(vs
に含まれない)しか通らないと仮定します.
(forall v, { in ps v, forall u, u \notin vs }) ->
ダイクストラ法はヒープを使って幅優先探索的なことをするアルゴリズムなので,最短距離が確定した頂点 u
は確定していない頂点 v
よりも早く見つかったのだぞという気持ちを込めて,u
への最短経路は v
へのいかなる経路より短いことを仮定します.
(forall u, u \notin vs -> { in vs, forall v p, last v p = s -> (cost_of_path c0 u (ps u) <= cost_of_path c0 v p)%O }) ->
ps v
は最短経路が確定した頂点を通って v
に向かう経路のうち最も短いものなので,それを仮定しているのですが,v
の最短距離が確定している場合は ps v
がグラフ全体を通して見ても最短経路であることも仮定しておきます.
(forall v p, { in p, forall u, u \notin vs } \/ v \notin vs -> last v p = s -> (cost_of_path c0 v (ps v) <= cost_of_path c0 v p)%O) ->
結論の部分は主定理と同様ですね.
(forall v, last v (dijkstra_rec c0 n vs ps v) = s) /\ (forall v, shortest_path c0 v (dijkstra_rec c0 n vs ps v)).
定理証明支援系を用いた形式検証では,帰納法が回るようなステートメントさえ決めてしまえば証明はすんなり行くのが良くあるパターンですが,この補題については一捻り必要でした.
問題となったのは,最短距離が確定していない頂点の中で最も短い経路 ps u
を持つ u
を選んだとき,最短経路の確定していない頂点へのどんな経路よりも ps u
は短いことの証明です.
have Hu : { in vs, forall v p, last v p = s -> (cost_of_path c0 u (ps u) <= cost_of_path c0 v p)%O }.
自然言語で何故これが成り立つかを解説しましょう.p
が最短経路の確定していない頂点を通ろうものなら,その時点で ps u
より長くなるので単調性から成り立ちます.また,ps
は最短経路の確定した頂点だけを通る経路の中で最短のものを覚えていて,最短距離が確定していない頂点の中で最も短い経路 ps u
を持つように頂点 u
を選んだので,p
が 最短経路の確定した頂点だけを通る場合も推移律から成り立ちます.
この途中で最短経路の確定していない頂点を通っていないことを言うのが大変で,ステートメントを一般化して帰納法を回してる最中にまたステートメントを一般化して帰納法を回すとか訳のわからないことをする羽目になりました.
have IHp : { in vs, forall v p p', { in p', forall u, u \notin vs } -> last (last v (rev p)) p' = s -> (cost_of_path c0 u (ps u) <= cost_of_path (cost_of_path c0 (last v (rev p)) p') v (rev p))%O }.
気持ちとしては,Hu
の p
を,最短距離の確定した頂点しか通らない前半の p’
と後半の p
に分けて,p
についての帰納法を回して p'
の占める範囲を増やしていき,実は最短距離の確定した頂点しか通らないと考えても良いことを言うイメージです.
最終的に,ダイクストラ法の正当性の証明は次のようになりました.
Lemma dijkstra_rec_correct c0 s : forall n vs, size vs = n -> uniq vs -> forall ps, (forall u, last u (ps u) = s) -> (forall v, { in ps v, forall u, u \notin vs }) -> (forall u, u \notin vs -> { in vs, forall v p, last v p = s -> (cost_of_path c0 u (ps u) <= cost_of_path c0 v p)%O }) -> (forall v p, { in p, forall u, u \notin vs } \/ v \notin vs -> last v p = s -> (cost_of_path c0 v (ps v) <= cost_of_path c0 v p)%O) -> (forall v, last v (dijkstra_rec c0 n vs ps v) = s) /\ (forall v, shortest_path c0 v (dijkstra_rec c0 n vs ps v)). Proof. elim => /= [ ? /size0nil -> ? | ? IHn [ | v vs ] //= /succn_inj ? ? ] ps Hs Hps Hsep Hlb. - split => // ? ?. rewrite Hs. exact /Hlb /or_intror. - set u := [arg min_(u < v in v :: vs) cost_of_path c0 u (ps u)]%O. have -> : (if v == u then vs else v :: rem u vs) = rem u (v :: vs) by []. have [ Hin Hlb' ] : u \in v :: vs /\ { in v :: vs, forall w, (cost_of_path c0 u (ps u) <= cost_of_path c0 w (ps w))%O }. { have /(@Order.TotalTheory.arg_minP _ _ _ v (fun u => u \in v :: vs) (fun u => cost_of_path c0 u (ps u))) : v \in v :: vs by rewrite in_cons eqxx. by inversion 1. } have IHp : { in v :: vs, forall x p p', { in p', forall u, u \notin v :: vs } -> last (last x (rev p)) p' = s -> (cost_of_path c0 u (ps u) <= cost_of_path (cost_of_path c0 (last x (rev p)) p') x (rev p))%O }. { move => x ?. elim => /= [ ? ? ? | y p IHp p' Hnotin ]. - apply /(@le_trans _ _ (cost_of_path c0 x (ps x))). + exact /Hlb'. + apply /Hlb; eauto. - rewrite rev_cons last_rcons => ?. case (boolP (y \in v :: vs)) => ?. + apply /(@le_trans _ _ (cost_of_path c0 y (ps y))). { exact /Hlb'. } apply /(@le_trans _ _ (cost_of_path c0 y p')). { by refine (Hlb _ _ (or_introl _) _). } exact /cost_of_path_increase. + rewrite cost_of_path_rcons. apply /(IHp (y :: p')) => // ?. by rewrite in_cons => /orP [ /eqP -> | /Hnotin ]. } have Hu : { in v :: vs, forall x p, last x p = s -> (cost_of_path c0 u (ps u) <= cost_of_path c0 x p)%O }. { move => ? /IHp IH p. move: (IH (rev p) [::]). rewrite /= revK. by apply. } apply /IHn => [ | | ? | w | w | w p ]. + by rewrite size_rem. + exact /rem_uniq. + by rewrite !(fun_if (last _)) /= !Hs if_same. + case (cost_of_path c0 w (ps w) <= adj w u (cost_of_path c0 u (ps u)))%O => ?; rewrite rem_filter // mem_filter negb_and. * move => /Hps ->. by rewrite orbT. * rewrite in_cons => /orP [ /= -> // | /Hps -> ]. by rewrite orbT. + rewrite rem_filter // mem_filter negb_and => /orP [ /negbNE /eqP -> x | ? x ]; rewrite mem_filter => /andP [ ? ? ] ? ?; rewrite (fun_if (cost_of_path _ _)) -minEle le_minl. * by rewrite Hu. * by rewrite Hsep. + rewrite rem_filter // mem_filter negb_and (fun_if (cost_of_path _ _)) cost_of_path_cons -minEle le_minl. case (boolP (w \in v :: vs)) => [ Hin' | Hnotin ? ? ]. * rewrite Hin' orbF => [ [ | /negbNE /eqP -> /(Hu _ Hin) -> // ] ]. { case: p => [ ? ? | x p Hnotin ]. - rewrite Hlb //. by left. - have /Hnotin : x \in x :: p by rewrite in_cons eqxx. rewrite mem_filter negb_and => /= /orP [ /negbNE /eqP -> | ? ] ?. + exact /orP /or_intror /adj_monotone /Hu. + apply /orP /or_introl /(@le_trans _ _ (cost_of_path c0 w (x :: ps x))). { apply /Hlb => //=. left => ?. by rewrite in_cons => /orP [ /eqP -> | /Hps ]. } apply /adj_monotone /Hlb => //. by right. } * apply /orP /or_introl /Hlb; eauto. Qed. Theorem dijkstra_correct c0 s : (forall u, last u (dijkstra c0 s u) = s) /\ (forall v, shortest_path c0 v (dijkstra c0 s v)). Proof. apply /dijkstra_rec_correct => [ | | u | v ? | ? | v [ /= ? <- | /= u p ] ]. - by rewrite size_rem ?cardE ?mem_enum. - apply /rem_uniq /enum_uniq. - rewrite !(fun_if (last _)) /=. by case (eqVneq u s). - case (eqVneq v s) => //. by rewrite inE (rem_filter _ (enum_uniq _)) mem_filter negb_and => /= ? ->. - rewrite (rem_filter _ (enum_uniq _)) mem_filter negb_and mem_enum orbF => /= /negbNE -> ? ? ? ?. exact /cost_of_path_increase. - by rewrite eqxx lexx. - rewrite (fun_if (cost_of_path _ _)) /= (rem_filter _ (enum_uniq _)) mem_filter negb_and mem_enum orbF /=. case (eqVneq v s) => /= [ -> ? ? | ? [ ] // Hnotin ]. + exact /(le_trans (adj_increase _ _ _)) /adj_monotone /cost_of_path_increase. + have /Hnotin : u \in u :: p by rewrite in_cons eqxx. rewrite mem_filter negb_and mem_enum orbF => /negbNE /eqP -> ?. exact /adj_monotone /cost_of_path_increase. Qed.
今回検証したダイクストラ法の実装は検証を楽にするために計算量を無視した面が強いのですが,それでも証明が長くなってしまいましたね.
まとめ
本記事では,頂点に辿り着くまでのコストに頂点から頂点に伸びる辺を通った際のコストを加えたものが関数として与えられた際に,が単調増加関数かつを満たせばダイクストラ法が最短経路を返すことを解説しました.
加えて,ナイーブな実装ではありますが,そのようなグラフに対してのダイクストラ法が最短距離を返すことを Coq で検証しました.コードの全文は GitHub に上げているので,手元で動かしたい方はそちらをどうぞ.
特殊なグラフについての最短距離を求める問題というと,どうしても E - Two Currencies のように頂点を拡張して状態の情報を載せてしまいがちですが,本記事で紹介した条件を満たせば辺のコストが可変であっても良いと頭の片隅に置いておくと,計算量を余計に増やさずにすむかもしれません.