2014-9 - Researchmap

THE UNIVERSITY OF TOKYO
数理情報工学演習第一Ｃ
プログラミング演習 (第9回 )
2014/06/09
DEPARTMENT OF MATHEMATICAL INFORMATICS
1
THE UNIVERSITY OF TOKYO
今日の内容:
 最短路問題（ダイクストラ法）
 課題：ダイクストラ法
 最終レポート：行列積の実行時間を調べる
2
THE UNIVERSITY OF TOKYO
最短路問題とダイクストラ (Dijkstra) 法
有向グラフ G = (V, E) を考える
V: 節点集合, 節点数 n
E: 枝集合, 枝数 m
グラフ G の各枝 (i,j)  E は長さ aij を持つ
ある節点 s  V から別の節点 t  V への路の中で，
最も長さの短いものを見つける問題を最短路問題という
15
2
50
30
20
10
1
5
80
3
4
3
15
THE UNIVERSITY OF TOKYO
下のネットワークで節点 1 を s, 節点 5 を t とすると，s から t へは 1245,
1235,1345,
12345,135 の5つの路がある
路の長さはそれぞれ 95, 85, 120, 110, 95 なので，
最短路は 1235 である
大規模なネットワークでは全ての路を数え上げる方法は実用的ではない
15
2
50
30
20
10
1
5
80
4
4
3
15
THE UNIVERSITY OF TOKYO
最適性の原理
節点 s から t への最短路 P が得られているとする
路 P に含まれる節点を1つ任意に選び, r とする
路 P は s から r までの部分 P1 と， r から t への
部分 P2 に分割できる．このとき，P1 は s から r への最短路で，P2 は r から t
への最短路となる
もし P1 より短い s から r への路 P1’が存在するなら，
P1’と P2 をつないだ路は P2 より短くなる
このような性質を最適性の原理と呼ぶ
P1
s
5
r
P1’
P2
t
THE UNIVERSITY OF TOKYO
ある節点 s からネットワークの他の全節点への最短路を求める問題を考える
ダイクストラ法は，枝の長さに関する非負条件
aij  0 ((i,j)  E) の仮定の下で，節点 s から各節点 i  V への最短路の長さ
の上限値 d(i) を更新していく反復アルゴリズム
計算の途中では，d(i) の値が s から i への真の最短路の長さに等しいことがわかっ
ている節点が存在する．そのような節点の集合を S で表す
S
6
は S の補集合 V  S を表す
THE UNIVERSITY OF TOKYO
ダイクストラ (Dijkstra) 法
(0) S := ,
とおく
S := V,
d(s) := 0, d(i) :=  (i  V  {s})
S = V なら計算終了．そうでないなら，d (v)
である節点 v  S を選ぶ
(2)
 min d (i) | i  S 
S : S  {v}, S : S  {v} とし，(v,j) E かつ j  S
であるような全ての枝 (v,j) に対して
d(j) > d(v) + avj ならば d(j) = d(v) + avj, p(j) := v
とする．ステップ (1) に戻る
p(i) は，s から i までの最短路において i の直前に
位置する節点を示すために用いられる
7
THE UNIVERSITY OF TOKYO
[初期化]
(0)
S   , S  {1,2,3,4,5}
50
1
d(1)=0
8
d(2)=
2
20
80
d(4)=
4
15
30
10
3
d(3)=
d(5)=
5
15
THE UNIVERSITY OF TOKYO
[反復1]
(1)
(2)
S  {1,2,3,4,5}
min d (i ) | i  S   0 より v = 1
S  {1}, S  {2,3,4,5}
d(2) =  > d(1) + a12 =50 より d(2) = 50, p(2) = 1
d(3) =  > d(1) + a13 =80 より d(3) = 80, p(3) = 1
50
1
d(1)=0
9
d(2)=50
2
20
80
d(4)=
4
15
30
10
3
d(3)=80
d(5)=
5
15
THE UNIVERSITY OF TOKYO
[反復2]
(1)
S  {2,3,4,5}
min d (i) | i  S   50 より v = 2
(2) S  {1,2}, S  {3,4,5}
d(3) = 80 > d(2) + a23 =70 より d(3) = 70, p(3) = 2
d(4) =  > d(2) + a24 =65 より d(4) = 65, p(4) = 2
50
1
d(1)=0
10
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
d(5)=
5
15
THE UNIVERSITY OF TOKYO
[反復3]
(1)
(2)
S  {3,4,5}
min d (i) | i  S   65 より v = 4
S  {1,2,4}, S  {3,5}
d(5) =  > d(4) + a45 =95 より d(5) = 95, p(5) = 4
50
1
d(1)=0
11
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=95
15
THE UNIVERSITY OF TOKYO
[反復4]
(1)
(2)
S  {3,5}
min d (i) | i  S   70 より v = 3
S  {1,2,3,4}, S  {5}
d(5) = 95 > d(3) + a35 =85 より d(5) = 85, p(5) = 3
50
1
d(1)=0
12
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=85
15
THE UNIVERSITY OF TOKYO
[反復5]
(1)
(2)
S  {5}
自動的に v = 5
S  {1,2,3,4,5}, S  
[反復6]
(1) S = V なので終了
50
1
d(1)=0
13
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=85
15
THE UNIVERSITY OF TOKYO
アルゴリズムが終了したときの d(i) は，節点 1 から
i への最短路の長さを与えている
枝の集合 {(p(i),i)} が各節点への最短路を表す
これを最短路木と呼ぶ
50
1
d(1)=0
14
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
d(5)=85
5
15
THE UNIVERSITY OF TOKYO
ダイクストラ法の実装

• d (v)  min d (i ) | i  S
d(v) をヒープに入れる
 である点 v を高速に求めるために
•
しかし，「d(j) > d(v) + avj ならば d(j) = d(v) + avj」を行うには，
ヒープ中の d(j) を更新し，ヒープを作り直す (heapify) 必要がある
•
d(j) は最小値ではないので，ヒープの配列のどこにあるか分からない
•
解決法
1. 点番号 j からヒープ中の位置を求める配列を追加する
2. ヒープ中の d(j) を更新しないようにアルゴリズムを変更する
15
THE UNIVERSITY OF TOKYO
解決法1
• ヒープに格納する値は (i, d(i)) で，1  i  n で i は相異なるとする
• ヒープの配列 A[1],…,A[n] で A[j] = (i, d(i)) のとき，
INV[i] = j とする (1  i  n )
• (i, d(i)) がヒープ内にないときは INV[i] = -1 などとする
• ヒープ内で要素を入れ替えたときは INV の値も更新する
16
THE UNIVERSITY OF TOKYO
解決法2
• d の値を更新するとき，古い値はそのままヒープに入れておき，
新しい値も別の要素として追加する
• d(j) = d(v) + avj と更新するとき，ヒープには (j, d(j)) が入っている
これの他に (j, d(v) + avj) もヒープに入れる
• ヒープ内の j には重複がある
• ヒープのサイズは枝数だけ必要 (解決法1では点数)
• ヒープから最小値を取り出すときに古い d の値が出ることがあるので，
そのときは無視する（後述）
17
THE UNIVERSITY OF TOKYO
初期状態
H = {(1, 0)}
ヒープから (1,0) を取り出す
50
1
d(1)=0
18
d(2)=
2
20
d(4)=
4
15
30
10
5
80
3
d(3)=
d(5)=
15
THE UNIVERSITY OF TOKYO
1から出ている枝で行ける点 2,3 について d を更新
(2,50), (3,80) をヒープに入れる
50
1
d(1)=0
19
d(2)=50
2
20
d(4)=
4
15
30
10
5
80
3
d(3)=80
d(5)=
15
THE UNIVERSITY OF TOKYO
H = {(2,50), (3,80)}
ヒープから (2,50) を取り出す
2から出ている枝で行ける点 3,4 について d を更新
(3,70), (4,65) をヒープに入れる
50
1
d(1)=0
20
d(2)=50
2
20
d(4)=65
4
15
30
10
5
80
3
d(3)=70
d(5)=
15
THE UNIVERSITY OF TOKYO
H = {(3,80), (3,70), (4,65)} (3が重複していることに注意)
ヒープから (4,65) を取り出す
4から出ている枝で行ける点 5 について d を更新
(5,95) をヒープに入れる
50
1
d(1)=0
21
d(2)=50
2
20
d(4)=65
4
15
30
10
5
80
3
d(3)=70
d(5)=
15
THE UNIVERSITY OF TOKYO
H = {(3,80), (3,70), (5,95) }
ヒープから (3,70) を取り出す
3から出ている枝で行ける点 4,5 について d を更新 (5のみ)
(5,85) をヒープに入れる
50
1
d(1)=0
22
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=85
15
THE UNIVERSITY OF TOKYO
H = {(3,80), (5,95), (5,85)}
ヒープから (3,80) を取り出す
しかし d(3) = 70 < 80 なので 3 への最短路はすでに求まっている
50
1
d(1)=0
23
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=85
15
THE UNIVERSITY OF TOKYO
H = {(5,95), (5,85)}
ヒープから (5,85) を取り出す
5 から出ている枝は無いので何もしない
50
1
d(1)=0
24
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=85
15
THE UNIVERSITY OF TOKYO
H = {(5,95)}
ヒープから (5,95) を取り出す
しかし d(5) = 85 < 95 なので 5 への最短路はすでに求まっている
50
1
d(1)=0
25
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=85
15
THE UNIVERSITY OF TOKYO
H = {}
ヒープが空なので終了
50
1
d(1)=0
26
d(2)=50
2
20
80
d(4)=65
4
15
30
10
3
d(3)=70
5 d(5)=85
15
THE UNIVERSITY OF TOKYO
課題
• ヒープを用いてダイクストラ法を実装する
• 解決法1,2のどちらでも良い
• 6/13(金) の正午までに提出
• 宛先：[email protected]
27
THE UNIVERSITY OF TOKYO
最終レポート
• 疎行列の2乗を計算し，実行時間を測定する
• 入力は Matrix Market のデータを用いる
http://math.nist.gov/MatrixMarket/matrices.html
– 形式を変換したものを用意しています
http://researchmap.jp/mub2jpgrp-1587/#_1587
• ファイルから行列を読み込み，その2乗を計算しファイルに出力する
– ./a.out (入力) (出力)
• 実行時間は time コマンドで測る
– time ./a.out bcspwr01.txt bcspwr01-2.txt
28
THE UNIVERSITY OF TOKYO
• 実行時間のグラフを作る
• 行列のデータ構造や乗算アルゴリズムは複数用いてそれらの違いを調べる
• それぞれのデータ構造，アルゴリズムの概略も説明する
• 締切：6/27(金) 正午
• 宛先：[email protected]
29
THE UNIVERSITY OF TOKYO
前回の補足
• 汎用的なヒープの実装として，2つ紹介した
– 方法1: ハッシュ関数のときと同じように，外部の比較関数を使う
– 方法2: マクロを使う
• 実装例は以下のページからダウンロードできます
– http://researchmap.jp/mu84xo6ez-1587/#_1587
30
THE UNIVERSITY OF TOKYO
方法1：外部の関数を用いる方法
• ヒープの構造体は heap.h で次のように定義する
typedef struct heapdata_ heapdata;
typedef struct heap_ {
int max_size;
int size;
heapdata *A; // void *A; としても良い
} HEAP;
• 外部関数として，次のものを使う
// A[i] > A[j] なら1を返す外部関数
int heap_greater(HEAP *H, int i, int j);
// A[i] と A[j] を入れ替える外部関数
void heap_swap(HEAP *H, int i, int j);
31
THE UNIVERSITY OF TOKYO
• 配列の中身は整数，実数，構造体など様々でサイズも違うので，
heap.c 内でデータの移動などは実現しにくい
• HEAPの構造体では，配列 A のポインタを格納しているだけで，
heap.c では配列の中身には触れないことにする
• heap_extract_max() や heap_insert() は heapdata を入出力にしていた
ので変更する
– extract_max は最小値を配列の最後に移動するだけ
– insert は，挿入する値を配列の最後に書き込んでから呼び出す
32
THE UNIVERSITY OF TOKYO
方法2: マクロを使う
• heap.h の最初に，データの宣言と比較方法を書く
typedef struct {
int x;
double y;
} heapdata;
#define GREATER(a, b) ((a).x > (b).x)
• それ以外の部分はデータの型に依存せずに書ける
– (ただし型を変えるときは再コンパイルが必要)
• 方法1,2 でも，1つのソース内で複数種類のヒープを用いることはできない
33
THE UNIVERSITY OF TOKYO
方法3: C++のテンプレート的手法
• ヒープに格納するデータごとに異なる構造体，プログラムを用意する
typedef struct {
int max_size;
int size;
int *A;
} heap_int;
heap_int *heap_build_int(int n, int *A, int max_size);
typedef struct {
int max_size;
int size;
string *A;
} heap_string;
heap_string *heap_build_string(int n, string *A, int max_size);
• 毎回書くのは大変なのでマクロを活用する
34
THE UNIVERSITY OF TOKYO
• ヒープを使いたいソースファイルに次のように書く
#define TYPE int
#define ___(s) s ## _int
#include "heap.h"
• heap.h には次のように書く
typedef struct {
int max_size;
int size;
TYPE *A;
} ___(HEAP);
• すると，heap.h に次のように書いたのと同じになる
typedef struct {
int max_size;
int size;
int *A;
} heap_int;
35
THE UNIVERSITY OF TOKYO
• ヒープの処理をする関数自体を heap.h 内に書く (プロトタイプ宣言では無く)
static TYPE ___(heap_extract_max)(___(HEAP) *H)
{
TYPE MAX, *A;
A = H->A;
if (H->size < 1) {
printf("heap_extract_max: ERROR∖n");
exit(1);
}
MAX = A[1];
A[1] = A[H->size];
H->size = H->size - 1;
___(heap_heapify)(H,1);
return MAX;
}
36
THE UNIVERSITY OF TOKYO

Download Report