最短経路問題に対するコストスケーリング法

概要

$n$ 頂点 $m$ 辺のグラフと(非負とは限らない)辺費用 $c_e$ が与えられる. この入力での最短経路問題は Bellman-Ford 法や SPFA で $O(nm)$ で計算可能であることがよく知られている. この記事では, 辺費用が $-C$ 以上の整数である場合に適用できる $O(\sqrt{n} m \log(C))$ のアルゴリズムをこの論文より紹介する. このアルゴリズムの一種は最小費用流問題に対するコストスケーリング法の最適化で, 具体的にはポテンシャルのみを変更することで scaling phase を飛ばすことが出来るか判定するために使われることがある. 個人的には, Dilworthの定理/Mirskyの定理の系を活用しているあたりが非常に面白いと思った.

準備

この節では, 後に必要となる chain と antichain を定義し, その大きさについての定理を示す.

$n$ 頂点の DAG (Directed Acyclic Graph, 有向非巡回グラフ) $G$ を考える. 頂点 $u$ から頂点 $v$ への(空でない)パスが存在するとき $u \succ v$, そうでないとき $u \not\succ v$ と書くことにする. 頂点部分集合 $S$ であって, 任意の相異なる2頂点 $u, v \in S$ が $u \succ v$ 又は $v \succ u$ を満たすものを chain (鎖) と呼ぶ. また, 頂点部分集合 $S$ であって, 任意の相異なる2頂点 $u, v \in S$ が $u \not\succ v$ かつ $v \not\succ u$ を満たすものを antichain (反鎖) と呼ぶ. 例えば, $G$ のパスの頂点集合は chain であり, 入次数 $0$ の頂点全体からなる集合は anti chain である.

chain と antichain の要素数の意味での最大値に関して, 次の定理が成り立つ:

定理1:

証明:

応用

問題設定と用語の準備

この辺りのことは ながたかなさんの記事 にも詳しく書かれているので, ポテンシャルや被約費用を扱ったことがない人は一読するとよいだろう.

$n$ 頂点 $m$ 辺の単純⁴有向グラフ $G = (V, E)$ と, 辺費用 $c(e) \in \Z$ が与えられる. 頂点から実数への関数 $\pi: V \to \R$ をポテンシャルと呼ぶ. ポテンシャル $\pi$ に対し, $c^\pi((u, v)) := c((u, v)) + \pi(u) - \pi(v)$ を 被約費用 (又は簡約費用, reduced cost) と呼ぶ.⁵ 任意の辺に対して被約費用が非負であるとき, そのポテンシャルはfeasible (実行可能) であると言う.

任意の $s$-$t$ パス $P$ と閉路 $C$ に対し, $$ \begin{aligned} c(P) &:= \sum_{e \in P} c(e), &c^\pi(P) &:= \sum_{e \in P} c^\pi(e), \\ c(C) &:= \sum_{e \in C} c(e), &c^\pi(C) &:= \sum_{e \in C} c^\pi(e) \end{aligned} $$ とすると, $c^\pi$ の定義から, $$ \begin{aligned} c^\pi(P) &= \sum_{(u, v) \in P} \left( c((u, v)) + \pi(u) - \pi(v) \right) \\ &= c(P) + \pi(s) - \pi(t), \\ \\ c^\pi(C) &= \sum_{(u, v) \in C} \left( c((u, v)) + \pi(u) - \pi(v) \right) \\ &= c(C) \end{aligned} $$ である. 従って, 閉路の被約費用はポテンシャルのとり方に依らない. また, パスの被約費用は, そのパス自体の費用, 始点と終点でのポテンシャルの値のみに依存する. よって, $s$-$t$ 最短経路は, その辺費用を $c$ としても $c^\pi$ としても変わらない.

ポテンシャル $\pi$ が feasible であるとき, その定義より $c^\pi(e) \ge 0$ であるから, 非負費用 $c^\pi$ による最短経路問題を解けばもとの費用 $c$ での最短経路問題と同じ解を得ることが出来る. 非負費用であるとき, 最短経路問題は Dijkstra 法で $O(m + n \log n)$ で解けることが知られているから, feasible なポテンシャルを高速に求めるアルゴリズムが欲しい.

上の議論では, feasible なポテンシャルの存在については言及していなかった. 例えば負閉路, すなわち $c(C) < 0$ となる閉路 $C$ があるとき, 任意のポテンシャル $\pi$ には, $c(C) = c^\pi(C) = \sum_{e \in C} c^\pi(e) < 0$ であるから $c^\pi(e) < 0$ となる辺が存在し, feasible でない. しかし, この時はそもそも最短経路が存在しない⁶. 一方, 負閉路が無ければ feasible なポテンシャルが存在する. 実際, $G$ に頂点 $s$ を加え, $s$ 以外の任意の頂点へ $s$ から費用 $0$ の辺を張ったグラフを考える. $s$ から $v$ への最短距離を $\pi(v)$ とすると, 任意の $G$ の辺 $e = (u, v)$ に対し, $\pi$ の最適性から $\pi(v) \le c(e) + \pi(u)$ であるが, これは $\pi$ が feasible なポテンシャルであることに他ならない.

以下の節での目的は, feasible なポテンシャル, すなわちポテンシャル $\pi: V \to \R$ であって, 任意の辺 $e = (u, v) \in E$ で $c^\pi(e) = c(e) + \pi(u) - \pi(v) \ge 0$ となるものが存在するか判定し, 存在する場合は実際に構成するアルゴリズムを作ることである.

緩和とスケーリング法

feasible なポテンシャルの定義は $\forall e \in E, c^\pi(e) \ge 0$ であった. $\epsilon > 0$ でこれを緩和した $\forall e \in E,\ c^\pi(e) > -\epsilon$ を満たすポテンシャル $\pi$ を $\epsilon$-feasible ($\epsilon$-実行可能) であるという. ⁷ また, $c^\pi(e) \le 0$ である辺を admissible (使用可能, 又は許容) であると言い, $c^\pi(e) \le -\epsilon$ である辺を $\epsilon$-improvable な辺 であると言う. ⁸ また, $\epsilon$-improvable な辺 $(u, v)$ が存在するような頂点 $v$ を $\epsilon$-improvable な頂点 であると言う.

入力の整数性と前節の feasible なポテンシャルの構成から, feasible なポテンシャルが存在することと feasible な整数ポテンシャルが存在することは同値であるから, 以下補足の節まで, ポテンシャルは整数ポテンシャルのみを考える. すると, $\epsilon = 1$ のとき, $\epsilon$-feasible なポテンシャル $\pi$ は $c^\pi((u, v)) = c((u, v)) + \pi(u) - \pi(v) > -1$ を満たすが, 左辺は整数であるから $c^\pi(e) \ge 0$, つまり feasible である.

さて, $C = \max(2, \max_e(-c_e))$ とする. 従って, $c_e \ge -C$ かつ $\log C > 0$ である⁹. $\epsilon = 2^{\floor{\log C} + 1}$ とすると, $c_e = c_e + 0 - 0 \ge -C > -\epsilon$ であるから, ゼロポテンシャル $\pi(\cdot) = 0$ は $\epsilon$-feasible である. コストスケーリング法は, この $\epsilon$ とゼロポテンシャルを初期値とし, より小さい $\epsilon$ に対する $\epsilon$-feasibility を満たすよう $\pi$ を変更していく, 次のようなアルゴリズムである.

1. $\epsilon = 2^{\floor{\log C} + 1}$, $\pi(\cdot) = 0$ とする.
2. 不変条件: $\pi$ は $\epsilon$-feasible で, $\epsilon$ の整数倍のみを値に持つポテンシャルである.
3. $\epsilon = 1$ であるならば終了し, そうでないならば $\epsilon := \epsilon / 2$ と更新する.
4. $\operatorname{Refine}(\epsilon, \pi)$ を呼び出す. この Refine は, $(2\epsilon)$-feasible なポテンシャルを $\epsilon$-feasible なポテンシャルに変形するか, その非存在を報告して(アルゴリズム全体を)終了するサブルーチンである.
5. 2.へ戻る.

Refine は $O(\log C)$ 回呼び出されるから, Refine の時間計算量が $O(f(n, m))$ であるならば, コストスケーリング法全体の時間計算量は $O(f(n, m) \log C)$ である. 以下, この Refine の実現を目標とする.

Admissible graph と Decycle

ポテンシャル $\pi$ を一つ固定する. このポテンシャルについて Admissible な辺 $E^\pi := \set{e \in E \setmid c^\pi(e) \le 0}$ のみを取り出したグラフ $G^\pi = (V, E^\pi)$ を admissible graph と呼ぶ. この $G^\pi$ の任意の閉路 $C$ のコストは $c(C) = c^\pi(C) = \sum_{e \in C} c^\pi(e) \le 0$ であり, $C$ が負の被約費用を持つ辺を含むならば $C$ は負閉路である. 従って, $G^\pi$ に両端点が同じ強連結成分に属する辺で負の被約費用を持つものが存在するかを判定し, 存在するならばアルゴリズムを終了してよい.

一方, そのような辺が存在しない場合を考える. 両端点が同じ強連結成分に含まれる辺全体の集合を $S$ とする. $G$ をもとにし, 強連結成分たちを縮約したグラフ $H$ を考えると¹⁰, $H$ 上のポテンシャル $\rho$ と $G$ 上の $S$ に含まれる辺の被約費用が $0$ であるようなポテンシャル $\rho'$ は, その整数性, feasibility や $\epsilon$-feasibility を保ったまま 1:1 対応し, $H$ の $\pi$ に対応するポテンシャルでの admissible graph は DAG になる. 従って, 必要ならば強連結成分分解をすることで, admissible graph は DAG であるとしてよい. この操作を decycle と呼ぶ. ¹¹

Cut-Relabel

さて, $(2\epsilon)$-feasible なポテンシャル $\pi$ と, (必要ならば decycleし) 付随する DAG である admissible graph $G^\pi$ を考える. 頂点集合 $S$ に対し, $G^\pi$ 上で $S$ から辿り着ける頂点全体を $S$ の閉包と呼び, $\overline{S}$ で表す. $S = \overline{S}$ であるとき, $S$ は閉集合であるという. 任意の頂点集合 $S$ に対し, $\overline{S}$ は閉集合である.

$S$ が閉集合であるとき, その要素に対応する $\pi$ を $\epsilon$ 減少させる操作を $\operatorname{CUT\_RELABEL}(S)$ で表す. この操作により, $c^\pi((u, v)) = c((u, v)) + \pi(u) - \pi(v)$ は $u \not \in S$ かつ $v \in S$ のとき $\epsilon$ 増加し, $u \in S$ かつ $v \not \in S$ のときに $\epsilon$ 減少する. しかし $S$ の選び方から, $u \in S$ かつ $v \not \in S$ であるような辺 $e = (u, v)$ は admissible でないから, この操作の後も $c^\pi(e) > 0 - \epsilon = - \epsilon$ を満たす. 従って, この操作は $\epsilon$-improvable な辺や頂点を増やさない. また, 更新前の $\pi$ が $\epsilon$ の整数倍であるならば更新後も $\epsilon$ の整数倍である.

どの二頂点間にも最後に $\epsilon$-improvable な辺を通るパスが存在しないような頂点部分集合をひとつ取り, $X$ とする. $e = (u, v)$ が $v \in X$ かつ $\epsilon$-improvable であるとすると, $X$ と閉包の定義から $u \not \in \overline{X}$ である. $\pi$ が $(2\epsilon)$-feasible であったことを思い出すと, この辺は $\operatorname{CUT\_RELABEL(\overline{X})}$ によって $c^\pi(e) > -2\epsilon + \epsilon = - \epsilon$ となり, $\epsilon$-improvable でなくなることがわかる. よって, このような頂点集合 $X$ に含まれる $\epsilon$-improvable な頂点は, $\operatorname{CUT\_RELABEL}(\overline{X})$ によって $\epsilon$-improvable でなくなる.

これにより, Refine の実現として, 次のアルゴリズムが考えられる.

1. 必要ならば Decycle する.
2. $\epsilon$-improvable な頂点が存在しない, つまり $\pi$ が $\epsilon$-feasible であれば, 終了する. そうでないとき, $\epsilon$-improvable な頂点をひとつ取り $x$ とする.
3. $\operatorname{CUT\_RELABEL(\overline{\set{x}})}$ を実行する
4. 1.へ戻る

上のアルゴリズムで, 2. を実行する度に $\epsilon$-improvable な頂点は減っていくから, 2. は高々 $n$ 回実行される. $\epsilon$-Improvable な頂点の発見や $\operatorname{CUT\_RELABEL(\overline{\set{x}})}$ の実行は $O(m)$ でできるから, 上の Refine は $O(nm)$ であり, コストスケーリング法全体の時間計算量は $O(nm \log C)$ となる.

Eliminate-Chain

ポテンシャル $\pi$ は $2 \epsilon$-feasible であり, その値は $\epsilon$ の整数倍であるとする. 任意の $x$ に対し, $x$ 以上の最小の $\epsilon$ の整数倍を $\epsceil{x}$ で表す. すなわち, $\epsceil{x}=\ceil{x/\epsilon}\epsilon$ とする. ¹² また, $c^\pi_\epsilon(e) = \epsceil{c^\pi(e)}$ とする. 例えば $c^\pi_\epsilon(P) = \sum_{e \in P} \epsceil{c^\pi(e)}$ である.

$P$ を $G^\pi$ のパスとする. パスにこの順で含まれる頂点 $u, v$ に対し, $P_v$ をこのパスの始点から $v$ までの部分パスとし, $P_{u,v}$ を $u$ から $v$ までの部分パスとする. 他のパスに対しても同様の記法を使う. さて, $$ \rho(v) = \begin{cases} c^\pi_\epsilon(P_v) &\text{if } v \in P, \\ 0 &\text{otherwise} \end{cases} $$ とすると, $P$ は $G^\pi$ のパスだったから, $\rho(v) \le 0$ であり, $\rho$ は $P$ に沿って単調非増加である. また, $c^\pi(e) > -2\epsilon$ であるから, $c^\pi_\epsilon(e) \ge -\epsilon$ であり, $c^\pi_\epsilon(P_v) \ge - n \epsilon$ である.

$G$ に新たな頂点 $s$ を加え, $s$ から他の全頂点に辺を加えたグラフを $H$ とし, $H$ 上の辺重み $c'$ を $$ c'(e) = \begin{cases} \rho(v) + n \epsilon & \text{if } e = (s, v), \\ \max\set{0, c^\pi_\epsilon(e)} & \text{otherwise} \end{cases} $$ で定めると, これは非負であるから, $s$ から他の全頂点への $c'$ を辺コストとする最短経路長が存在する. これを $d(v)$ とする. 更に, $\pi'(v) = \pi(v) + d(v) - n \epsilon$ とする. この $\pi'$ によって $\pi$ を更新するのが Eliminate-Chain である.

天下り的に定義したが, これは $P$ に含まれる $\epsilon$-improvable な辺の終点について順に $\operatorname{CUT\_RELABEL}(\overline{\set{x}})$ を行った結果を求めていると考えてよい. 実際 $c^\pi_\epsilon(P_v)$ は $P$ の始点から順に $\epsilon$-improvable な辺を通るごとに $-\epsilon$ され, $c'((s,v))$ はこれを表している. また, $u$ が $\max\set{0, c^\pi_\epsilon((u,v)) / \epsilon}$ 回 relabel されると辺 $(u, v)$ は admissible になり, それ以降 $u$ を relabel する度に $v$ も relabel されることを $c'(u, v)$ は表している.

さて, $c'$ は $\epsilon$ の整数倍であり, 最短経路長は $d(v) \le c'((s, v)) = \rho(v) + n \epsilon \le n \epsilon$ を満たすから, Dial の方法¹³による Dijkstra 法を実行することにより $O(m)$ で計算可能である. また, $\pi$, $c'$ は $\epsilon$ の整数倍であるから, $d$, 従って $\pi'$ も $\epsilon$ の整数倍である.

次に, $\pi'$ は $(2\epsilon)$-feasible であり, 新たに $\epsilon$-improvable な辺が作られないことを示す. 任意の辺 $e = (u, v)$ を取る. $d$ の最適性から $c'(e) + d(u) - d(v) \ge 0$ であるから, $$ \begin{aligned} c^{\pi'}(e) &= c(e) + \pi'(u) - \pi'(v) \\ &= c(e) + (\pi(u) + d(u) - n \epsilon) - (\pi(v) + d(v) - n \epsilon) \\ &= c^\pi(e) + d(u) - d(v) \\ &= \left( c'(e) + d(u) - d(v) \right) + \left(c^\pi(e) - c'(e) \right) \\ &\ge c^\pi(e) - c'(e) \\ &= c^\pi(e) - \max\set{0, c^\pi_\epsilon(e)} \end{aligned} $$ となる. $e$ がもともと $\epsilon$-improvable であるとき $c^\pi(e) - \max\set{0, c^\pi_\epsilon(e)} = c^\pi(e) - 0 > - 2 \epsilon$ であり, $e$ がもともと $\epsilon$-improvable でなければ $c^\pi(e) - \max\set{0, c^\pi_\epsilon(e)} = c^\pi(e) - \epsceil{c^\pi(e)} > - \epsilon$ であるから, いずれにせよ $c^{\pi'}(e) > -2 \epsilon$ であり, $\pi'$ は $\pi$ から $\epsilon$-improvable な辺を増やさないことがわかった.

補題1:

定理2:

注意:

さて, 上の定理から, $\epsilon$-improvable な辺一つからなるパスを $P$ として Eliminate-Chain を繰返すことにより, $O(nm)$ で Refine を実現出来ることがわかる.

より高速な Refine の実現

Decycle した後の $G^\pi$ に対し, その $\epsilon$-improvable な頂点全体を $U$ とし, $k = |U|$ とする. 頂点集合を $U$ とし, $u, v \in U$ に対し, $u$ から $v$ へのパスで最後に $\epsilon$-improvable な辺を通るものが存在するときに辺を張ったグラフを $H$ とする.

$G^\pi$ が DAG であったから $H$ も DAG であり, 定理1 から $H$ の chain 又は antichain で要素数が $\sqrt{k}$ 以上のものが存在する. $H$ の chain $S \subset U$ は, $S$ を含む $G$ のパス $P$ で, $S$ へ入る辺が $\epsilon$-improvable なものへ展開することができる. 一方 $H$ の anti chain $T \subset V$ は, どの二頂点間にも最後に $\epsilon$-improvable な辺を通るパスが存在しない. 従って, $H$ の chain で要素数が $\sqrt{k}$ 以上のものがあればそれに Eliminate-Chain, そうでなければ $H$ の antichain で要素数が $\sqrt{k}$ 以上のものに Cut-Relabel を適用することにより, $O(m)$ で¹⁴少なくとも $\sqrt{k}$ 個の $\epsilon$-improvable な頂点を $\epsilon$-improvable でなくすことが出来る.

これを $i$ 回繰り返した後の $\epsilon$-improvable な頂点の数を $k_i$ とする. $k_i \ge k_0/2$ である間 $k_i - k_{i+1} \ge \sqrt{k_0 / 2}$ であるから, このような $i$ は高々 $(k_0 / 2) / \sqrt{k_0 / 2} = \sqrt{k_0 / 2}$ 回しかない. 同様にして $i$ 回目から $O(\sqrt{k_i/2})$ 回後には少なくとも半減するから, 合計の実行回数は $$ O\left(\sum_{i=0}^{\infty} \sqrt{k_0 / 2^i}\right) = O(\sqrt{n}) $$ である.

以上より, $O(\sqrt{n} m \log C)$ のアルゴリズムが得られた.

補足1: feasibility, improvable, admissible などを別の定義にした場合

$\epsilon$-feasibility を $c^\pi(e) \ge -\epsilon$ で定義することも出来る. 同時に, admissible edge は $c^\pi(e) < 0$, $\epsilon$-improvable edge は $c^\pi(e) < -\epsilon$ で定める. すると, admissible edge から $=0$ が抜けたことにより, admissible graph は (decycle せずとも) DAG になる.

一方で, $1$-feasible な整数ポテンシャルであっても feasible とは限らなくなる. しかし, $\pi$ が $\epsilon < \frac{1}{n}$ な $\epsilon$ に対する $\epsilon$-feasible なポテンシャルであるとき, $G$ の任意の閉路 $C$ のコストは $c(C) = c^\pi(C) = \sum_{e \in C} c^\pi(e) > n \cdot \frac{-1}{n} = -1$ であるから, 入力の整数性より非負である. 従って, feasible なポテンシャルが存在する.

以下, $\epsilon < \frac{1}{n}$, $\pi$ は $\epsilon$-feasible なポテンシャルであるとし, ここから feasible なポテンシャルを得る方法を考える.

$c'(e) = c^\pi(e) + \epsilon$ とすると, $\pi$ の $\epsilon$-feasibility から $c'(e) \ge 0$ である. また, 任意のパス $P$ に対し, $0 \le c'(P) - c^\pi(P) < n \epsilon < 1$ であるから, $c'$ での最短経路は $c^\pi$ での最短経路に一致する. 従って, 超頂点 $s$ と, $s$ から各頂点へ $c'((s, v)) = 0$ の辺 $(s, v)$ を加えたグラフ上で $s$ から各頂点への $c'$-最短経路木を求め, この最短経路木から $s$ を除いた森の上での被約費用が $0$ となるようポテンシャルを定めれば, これが feasible なポテンシャルになる.

従って, この方針では feasible なポテンシャルを合計 $O(\sqrt{n} m \log(n C))$ で計算することが出来る. 実際には, 入力のコストを $(n+1)$ 倍し, 常に整数に保つ実装にするとよい.

補足2: 実用上の話について

上では Refine の $O(nm)$ の実装を二つと $O(\sqrt{n} m)$ の実装を一つ紹介したが, 論文中には他にもいくつかの方法が与えられている. おそらく実用的には, REFINE-P の方が速いのではないかと思う. また, decycle, 特に縮約やそれを元に戻す操作をしなくてよいぶん, $O(\sqrt{n} m \log(nC))$ の方が実測で速くなる可能性も高いだろう.

参考文献

• Andrew V. Goldberg, 1993, "Scaling Algorithms for the Shortest Paths Problem", Proceedings of the fourth annual ACM-SIAM symposium on Discrete algorithms. Society for Industrial and Applied Mathematics, USA, 222-231.