PythonでAizu Online Judgeを解く/GRL_1_C - All Pairs Shortest Path
概要
重み付きの有向グラフ \(G(V, E)\) が与えられます。任意の2頂点について最短経路の距離を求めてください。また、グラフが負閉路を持つ場合はそれを判定しましょう。
制約
- \(1 \leq |V| \leq 100\)
- \(0 \leq |E| \leq 9900\)
- \(-2 \times 10^7 \leq d_i \leq 2 \times 10^7\)
- グラフ \(G\) は平行辺を持たない
- グラフ \(G\) は自己ループを持たない
入力
%python def input(f): global v, e, s, t, d v, e = map(int, f.readline().split()) s = [] t = [] d = [] for _ in range(e): s_, t_, d_ = map(int, f.readline().split()) s.append(s_) t.append(t_) d.append(d_)
サンプル入力
%sh cat << EOF > /tmp/input1 4 6 0 1 1 0 2 5 1 2 2 1 3 4 2 3 1 3 2 7 EOF cat << EOF > /tmp/input2 4 6 0 1 1 0 2 -5 1 2 2 1 3 4 2 3 1 3 2 7 EOF cat << EOF > /tmp/input3 4 6 0 1 1 0 2 5 1 2 2 1 3 4 2 3 1 3 2 -7 EOF
Python: 入力をnetworkxで可視化してみる
%sh apt update && apt install -y python3-pip pip3 install networkx matplotlib
WARNING: apt does not have a stable CLI interface. Use with caution in scripts. Hit:1 http://cran.rstudio.com/bin/linux/ubuntu xenial/ InRelease Hit:2 http://security.ubuntu.com/ubuntu xenial-security InRelease Hit:3 http://archive.ubuntu.com/ubuntu xenial InRelease Hit:4 http://archive.ubuntu.com/ubuntu xenial-updates InRelease Hit:5 http://archive.ubuntu.com/ubuntu xenial-backports InRelease Reading package lists... Building dependency tree... Reading state information... 69 packages can be upgraded. Run 'apt list --upgradable' to see them. WARNING: apt does not have a stable CLI interface. Use with caution in scripts. Reading package lists... Building dependency tree... Reading state information... python3-pip is already the newest version (8.1.1-2ubuntu0.4). 0 upgraded, 0 newly installed, 0 to remove and 69 not upgraded. Requirement already satisfied (use --upgrade to upgrade): networkx in /usr/local/lib/python3.5/dist-packages Collecting matplotlib Downloading https://files.pythonhosted.org/packages/7b/ca/8b55a66b7ce426329ab16419a7eee4eb35b5a3fbe0d002434b339a4a7b09/matplotlib-3.0.0-cp35-cp35m-manylinux1_x86_64.whl (12.8MB) Requirement already satisfied (use --upgrade to upgrade): decorator>=4.3.0 in /usr/local/lib/python3.5/dist-packages (from networkx) Collecting kiwisolver>=1.0.1 (from matplotlib) Downloading https://files.pythonhosted.org/packages/7e/31/d6fedd4fb2c94755cd101191e581af30e1650ccce7a35bddb7930fed6574/kiwisolver-1.0.1-cp35-cp35m-manylinux1_x86_64.whl (949kB) Collecting numpy>=1.10.0 (from matplotlib) Downloading https://files.pythonhosted.org/packages/75/22/355e68c80802d6f488223788fbda75c1daab83c3ef609153676c1f17be5f/numpy-1.15.2-cp35-cp35m-manylinux1_x86_64.whl (13.8MB) Collecting cycler>=0.10 (from matplotlib) Using cached https://files.pythonhosted.org/packages/f7/d2/e07d3ebb2bd7af696440ce7e754c59dd546ffe1bbe732c8ab68b9c834e61/cycler-0.10.0-py2.py3-none-any.whl Collecting pyparsing!=2.0.4,!=2.1.2,!=2.1.6,>=2.0.1 (from matplotlib) Downloading https://files.pythonhosted.org/packages/2b/4a/f06b45ab9690d4c37641ec776f7ad691974f4cf6943a73267475b05cbfca/pyparsing-2.2.2-py2.py3-none-any.whl (57kB) Collecting python-dateutil>=2.1 (from matplotlib) Using cached https://files.pythonhosted.org/packages/cf/f5/af2b09c957ace60dcfac112b669c45c8c97e32f94aa8b56da4c6d1682825/python_dateutil-2.7.3-py2.py3-none-any.whl Requirement already satisfied (use --upgrade to upgrade): setuptools in /usr/lib/python3/dist-packages (from kiwisolver>=1.0.1->matplotlib) Collecting six (from cycler>=0.10->matplotlib) Downloading https://files.pythonhosted.org/packages/67/4b/141a581104b1f6397bfa78ac9d43d8ad29a7ca43ea90a2d863fe3056e86a/six-1.11.0-py2.py3-none-any.whl Installing collected packages: kiwisolver, numpy, six, cycler, pyparsing, python-dateutil, matplotlib Successfully installed cycler-0.10.0 kiwisolver-1.0.1 matplotlib-3.0.0 numpy-1.15.2 pyparsing-2.2.2 python-dateutil-2.7.3 six-1.11.0 You are using pip version 8.1.1, however version 18.1 is available. You should consider upgrading via the 'pip install --upgrade pip' command.
%python import os if os.environ.get('ZEPPELIN_HOME'): import networkx as nx import matplotlib.pyplot as plt def visualize(): G = nx.MultiDiGraph() for i in range(v): G.add_node(i, node_color='b') for i in range(e): G.add_edge(s[i], t[i], t=d[i]) plt.figure() plt.style.use('seaborn') plt.axis("off") pos = nx.spring_layout(G) nodes = nx.draw_networkx_nodes(G, pos, node_color='w', linewidths=2) nodes.set_edgecolor('black') nx.draw_networkx_labels(G, pos) nx.draw_networkx_edges(G, pos) edge_labels = dict([((a, b), c['t']) for a, b, c in G.edges(data=True)]) nx.draw_networkx_edge_labels(G, pos, edge_labels=edge_labels, label_pos=0.25, font_size=16) plt.show() z.show(plt)
%python if os.environ.get('ZEPPELIN_HOME'): with open('/tmp/input1') as f: print('Example #1') input(f) visualize()
Example #1
%python if os.environ.get('ZEPPELIN_HOME'): with open('/tmp/input2') as f: print('Example #2') input(f) visualize()
Example #2
%python if os.environ.get('ZEPPELIN_HOME'): with open('/tmp/input3') as f: print('Example #3') input(f) visualize()
Example #3
頂点2と3の間を往復することでコストを下げ続けることができます。
考え方 - Warshall–Floyd Algorithm
100頂点くらいの規模感なグラフ上ですべての頂点のペアについて最短経路問題を解くとなると Warshall–Floyd 法を使え使えと言われているような気がします。このアルゴリズムでは2つの頂点について別のある頂点を経由したときに掛かるコストと現在の最小コストを比較して、別の頂点を経由したほうがコストが低ければ最小コストを更新していきます。
何となくで運用できてしまうアルゴリズムだと思うのですが、仕組みとしては動的計画法に基づいたものです。
参考: “ALGORITHM NOTE Warshall Floyd’s algorithm ワーシャル フロイド”
http://algorithms.blog55.fc2.com/blog-entry-52.html
$$
A^k[i, j] = 頂点 \{1, 2, 3, \dots, k\} のみを経由したときの頂点iと頂点jの最小コスト
$$
$$
P^k[i, j] = そのような最短経路の一つ
$$
とすると \(k=0\) のとき \(A^0[i, j]\) が最小コストであるのは明らかです。また \(k>0\) のとき、\(A^k\) について \(P^k[i, j]\) が頂点 k を経由しない場合は \(P^{k-1}[i, j]\) と等しくなり、頂点kを経由する場合は \(P^k[i, k]\) と \(P^k[k, j]\) の2つの部分路に分けることができてそれぞれ \(P^{k-1}[i, k]\) と \(P^{k-1}[k, j]\) に等しいものと考えることができます。したがって次のような漸化式を立てることができます。
$$
A^k[i,j] = min(A^{k-1}[i, j], A^{k-1}[i,k] + A^{k-1}[k,j])
$$
また、\(A^k[i, i] = 0\) より
$$
A^k[i,k] = min(A^{k-1}[i, k], A^{k-1}[i,k] + A^{k-1}[k,k]) = A^{k-1}[i,k]
$$
$$
A^k[k,j] = min(A^{k-1}[k, j], A^{k-1}[k,k] + A^{k-1}[k,j]) = A^{k-1}[k,j]
$$
となるので、実装においては各ステップで2次元配列の値を使い回すことができます。
実装
%python import sys def init(): global inf inf = sys.maxsize def generate_graph(): global g g = [[inf for _ in range(v)] for _ in range(v)] for i in range(v): g[i][i] = 0 for i in range(e): g[s[i]][t[i]] = d[i] def solve(): init() generate_graph() for k in range(v): for i in range(v): for j in range(v): if g[i][k] == inf or g[k][j] == inf: continue g[i][j] = min([g[i][j], g[i][k] + g[k][j]]) for i in range(v): if g[i][i] != 0: print('NEGATIVE CYCLE') return for i in range(v): for j in range(v): if j > 0: print(' ', end='') if g[i][j] == inf: print('INF', end='') else: print(g[i][j], end='') print('')
Example Output #1
%python if os.environ.get('ZEPPELIN_HOME'): os.environ['INPUT'] = '/tmp/input1' with open(os.environ.get('INPUT') or '/dev/stdin') as f: input(f) solve()
0 1 3 4 INF 0 2 3 INF INF 0 1 INF INF 7 0
Example Output #2
%python if os.environ.get('ZEPPELIN_HOME'): os.environ['INPUT'] = '/tmp/input2' with open(os.environ.get('INPUT') or '/dev/stdin') as f: input(f) solve()
0 1 -5 -4 INF 0 2 3 INF INF 0 1 INF INF 7 0
Example Output #3
%python if os.environ.get('ZEPPELIN_HOME'): os.environ['INPUT'] = '/tmp/input3' with open(os.environ.get('INPUT') or '/dev/stdin') as f: input(f) solve()
NEGATIVE CYCLE
%sh apt install -y jq ~/.aoj/bin/login curl -s http://localhost:8080/api/notebook/2DTWUY27U | \ jq -c -r '.body.paragraphs[].text | select(.|test("^%python"))' | \ sed 's/^%python//' | \ ~/.aoj/bin/submit GRL_1_C Python3 /dev/stdin exit 1
WARNING: apt does not have a stable CLI interface. Use with caution in scripts. Reading package lists... Building dependency tree... Reading state information... jq is already the newest version (1.5+dfsg-1ubuntu0.1). 0 upgraded, 0 newly installed, 0 to remove and 68 not upgraded. LOGIN_OK {"token": "614d410f-30bd-4996-b642-ada53856b0f8"}
計算量的には重いけれど実装はとてもシンプルで素早く問題を解決することができます。Warshall–Floyd 法はトレードオフという概念を初めて自分に突きつけてくれた思い入れのあるアルゴリズムの一つです。
%md