[PPT] in Japanese - 数理・計算科学専攻
Download
Report
Transcript [PPT] in Japanese - 数理・計算科学専攻
半正定値計画問題に対する
行列補完理論の高速実装
東京工業大学 数理・計算科学専攻
山下 真
この研究は笹川研究助成の助成を受けています
2012/12/15
SOTA:「最適化の理論と応用」研究部会
1
行列補完の簡単な例
半正定値計画問題(SemiDefinite Program, SDP)の例
Notation
内積計算に必要なのは青の要素のみ
これだけで計算できれば高速化できる
2012/12/15
SOTA:「最適化の理論と応用」研究部会
半正定値制約がある
2
組合せ最適化のときにあらわれる疎構造
格子ネットワーク(p,q)上の最大クリーク問題
内積計算に必要な要素のみを青で表示
ほとんどがゼロ要素⇒こういった特徴を利用して高速化
半正定値条件を取り入れる⇒行列補完
2012/12/15
SOTA:「最適化の理論と応用」研究部会
3
行列補完の簡単な例の続き
黒なしで計算
2012/12/15
SOTA:「最適化の理論と応用」研究部会
黒を適切に補完
4
今日の話の方向性
大規模な半正定値計画問題に対して
行列補完を利用して高速に求解
行列補完自体は SDPA-C に実装済み
今回のメインは SDPA-C の改良にあたる
話の内容の多くは、行列補完について
2012/12/15
SOTA:「最適化の理論と応用」研究部会
5
インデックス
1. 半正定値計画問題の標準形と主双対内点法
2. 行列補完理論と主双対内点法
3. 実装の高速化
1. 行列補完計算の改善
2. Schur complement マルチスレッド化
3. 主探索方向マルチスレッド化
4. 数値実験結果
2012/12/15
SOTA:「最適化の理論と応用」研究部会
6
半正定値計画問題の主な応用
制御理論
組合せ最適化問題の緩和
ロバスト最適化
量子化学における電子構造計算
量子計算
(物理状態は半正定値行列で表現できるため)
2012/12/15
SOTA:「最適化の理論と応用」研究部会
7
半正定値計画問題の標準形
(SemiDefinite Programs, SDP)
主問題と双対問題で定義
行列補完理論で本質的なのは双対側
2012/12/15
SOTA:「最適化の理論と応用」研究部会
8
主双対内点法の特徴
主問題と双対問題を同時に求解
多項式時間アルゴリズム
多くのソフトウェア
SDPA, CSDP, SeDuMi, SDPT3, …
主問題の変数が密行列
Schur 補完行列がボトルネック
2012/12/15
SOTA:「最適化の理論と応用」研究部会
9
主双対内点法の枠組み
実行可能領域
探索方向
ステップ長
中心パス
最適解
2012/12/15
SOTA:「最適化の理論と応用」研究部会
10
中心パスの点と探索方向計算
最適性条件
中心パスの点
探索方向
2012/12/15
の計算はニュートン法に基づく
SOTA:「最適化の理論と応用」研究部会
11
探索方向計算の続き
変数を左に集めて
変形していくと
まとめると
2012/12/15
SOTA:「最適化の理論と応用」研究部会
12
行列補完理論で
計算式が変わる部分
Schur 補完方程式
主問題側変数の計算
行列補完理論では密行列を避けて疎行列(下三角)で計算する
2012/12/15
SOTA:「最適化の理論と応用」研究部会
13
インデックス
1. 半正定値計画問題の標準形と主双対内点法
2. 行列補完理論と主双対内点法
3. 実装の高速化
1. 行列補完計算の改善
2. Schur complement マルチスレッド化
3. 主探索方向マルチスレッド化
4. 数値実験結果
2012/12/15
SOTA:「最適化の理論と応用」研究部会
14
行列補完理論と主双対内点法
変数行列の疎構造に基づいてコレスキー分解
Schur 補完行列を効率的に計算できる
詳しくは、以下の論文
Exploiting sparsity in semidefinite programming via matrix
completion I: general framework
Fukuda, Kojima, Murota, Nakata
SIAM J. Optim, 2000
Exploiting sparsity in semidefinite programming via matrix
completion II: implementation and numerical results
Nakata, Fujisawa, Fukuda, Kojima, Murota
Math. Program. B, 2003
SOTA:「最適化の理論と応用」研究部会
2012/12/15
15
行列補完理論のポイント
どれだけの要素が必要か⇒ Aggregate Sparsity Pattern
どう行列を分解するか
⇒ Chordal Graph & 極大クリーク
どう補完するか
⇒
2012/12/15
SOTA:「最適化の理論と応用」研究部会
16
Aggregate Sparsity Pattern
双対側の変数行列の非ゼロパターン
例
1
2
4
3
グラフ表現
5
7
6
2012/12/15
SOTA:「最適化の理論と応用」研究部会
17
Chordal Graph
グラフが Chordal であるとは、4以上のサイクルに弦(chord)があること
1
2
length4
4
3
length5
5
6
7
1
2
4
3
Chordal 化
5
7
6
ここで、極大 Clique は
(Clique は、すべての頂点間に枝がある部分のこと)
Point :: 極大 Clique に沿って、行列が分解される
2012/12/15
SOTA:「最適化の理論と応用」研究部会
18
行列の分解
1
2
4
3
5
7
6
青はAggregate, 赤は Chordal で追加
Grone et al. 1984
2012/12/15
SOTA:「最適化の理論と応用」研究部会
19
Chordal Graph と Cholesky分解
1
2
4
3
5
1
2
4
3
6
5
7
7
6
2012/12/15
SOTA:「最適化の理論と応用」研究部会
20
行列の分解
1
2
4
3
5
7
6
青はAggregate, 赤は Chordal で追加
Grone et al. 1984
黒の要素をどうする？
2012/12/15
SOTA:「最適化の理論と応用」研究部会
21
行列補完の簡単な例
1
正則転置
2012/12/15
正定値
SOTA:「最適化の理論と応用」研究部会
2
3
正則
22
行列補完のステップ１
1
2
4
3
5
7
6
2012/12/15
SOTA:「最適化の理論と応用」研究部会
23
行列補完のステップ２
下三角行列
これらを使うと、以下のように補完される
以下を満たすことを計算で確認できる
2012/12/15
SOTA:「最適化の理論と応用」研究部会
24
行列補完の性質
補完行列は密だが逆行列が疎構造を持つ
.
2012/12/15
SOTA:「最適化の理論と応用」研究部会
25
Cholesky分解に似た分解
Point::
2012/12/15
SOTA:「最適化の理論と応用」研究部会
26
双対側は通常のCholesky分解
Sparse Cholesky 分解を適用
2012/12/15
SOTA:「最適化の理論と応用」研究部会
27
行列補完と主双対内点法
主双対内点法のボトルネック
Schur 補完行列
主問題側の探索方向
密行列を避けて疎行列で計算する
2012/12/15
SOTA:「最適化の理論と応用」研究部会
28
Schur 補完行列の計算式の変更
2012/12/15
SOTA:「最適化の理論と応用」研究部会
29
主問題側探索方向の計算
各列ごとに計算
2012/12/15
SOTA:「最適化の理論と応用」研究部会
30
SDPA-C
行列補完理論を実装したSDPソルバー
SDPA(SemiDefinie Programming Algorithm) – Completion version
SDPA
従来の主双対内点法を実装
http://sdpa.sourceforge.net/ で公開
Debian/Ubuntu/Linux Mint パッケージ
$ sudo apt-get install sdpa libsdpa-dev sdpam
ソースコードは git でバージョン管理
４年分の全部の記録を管理
（複数サーバでの数値実験に便利）
2012/12/15
SOTA:「最適化の理論と応用」研究部会
SDPA と SDPA-Cで
ソースを比較すると
１万９千行のうち
１万行ぐらい違う
31
インデックス
1. 半正定値計画問題の標準形と主双対内点法
2. 行列補完理論と主双対内点法
3. 実装の高速化
1. 行列補完計算の改善
2. Schur complement マルチスレッド化
3. 主探索方向マルチスレッド化
4. 数値実験結果
2012/12/15
SOTA:「最適化の理論と応用」研究部会
32
SDPA-CとSDPAの比較
SDPA-C
SDPA-C6.2.1
6.2.1
SDPA
SDPA7.3.8
7.3.8
SDPA-C
SDPA-C7.3.8
7.3.8
行列補完理論
行列補完理論
SparseCholesky
Cholesky分解
分解
Sparse
○
○
自作ルーチン
自作ルーチン
×
×
MUMPS
MUMPS
Schurマルチスレッド
マルチスレッド
Schur
dXマルチスレッド
dXマルチスレッド
×
×
×
×
○
○
×
×
○
○
CHOLMOD
CHOLMOD
&MUMPS
&MUMPS
○
○
CallableLibrary
Library
Callable
Matlab対応
対応
Matlab
×
×
×
×
○
○
○
○
2012/12/15
SOTA:「最適化の理論と応用」研究部会
○
○
○
○
○
○
33
行列分解計算式の検討
.
.
Point
2012/12/15
SOTA:「最適化の理論と応用」研究部会
34
2012/12/15
SOTA:「最適化の理論と応用」研究部会
35
2012/12/15
SOTA:「最適化の理論と応用」研究部会
36
これと同じことを極大クリークの情報で
帰納法を行えば一般的にできる
2012/12/15
クリークの情報にアクセスできる
数値計算ライブラリ
SOTA:「最適化の理論と応用」研究部会
37
Sparse Cholesky 分解のライブラリ
CHOLMOD
アルゴリズムは supernodal 分解
行列演算のマルチスレッドで性能が伸びにくい
行列補完理論に採用
MUMPS
a MUltifrontal Massively Parallel sparse direct Solver
アルゴリズムは multiple frontal method
行列演算のマルチスレッドで高速化される
Schur complement 方程式の求解に採用
2012/12/15
SOTA:「最適化の理論と応用」研究部会
38
Supernode 分解
.
Sparse Cholesky 分解を効率的に行える
supernode から
クリークの情報が得られる
2012/12/15
SOTA:「最適化の理論と応用」研究部会
39
内点法の前処理
2012/12/15
SOTA:「最適化の理論と応用」研究部会
40
数式改善による高速化
格子上の最大クリーク問題
クリーク数 438, 平均29.89, 最大59, 合計 13090
SDPA-C6.2.0
SDPA-C7.3.8
Schur 要素計算
4205.70秒
2094.03秒
Schur Cholesky
185.36秒
161.03秒
dX の計算
316.00秒
242.51秒
22.22秒
17.93秒
その他
合計時間
4729.28秒
2515.50秒 1.88倍高速化
時間は短縮されたが、まだまだ遅い
マルチスレッドによる並列計算で高速化
2012/12/15
SOTA:「最適化の理論と応用」研究部会
41
Schur 補完行列のマルチスレッド化
各列の計算は独立している
スレッド割り当て⇒終わったスレッドが次の列へ
4
3
2
1
2012/12/15
2
3
4
1
このままだとスレッドが衝突する
SOTA:「最適化の理論と応用」研究部会
42
スレッドの衝突
行列積担当の BLAS もマルチスレッド
そのままマルチスレッド化すると（４コアのとき）
スレッド1
Schur 補完行列
スレッドが増えすぎて
スレッドの切り替えに
時間を取られてしまう
スレッド2
スレッド3
スレッド4
列計算スレッド
2012/12/15
BLASスレッド
SOTA:「最適化の理論と応用」研究部会
43
スレッドの衝突回避
Schur 補完行列計算直前にBLASのスレッドを無効化
Schur 補完行列計算終了後にBLASのスレッドを戻す
Schur 補完行列
2012/12/15
スレッド1
スレッド1
スレッド2
スレッド2
スレッド3
スレッド3
スレッド4
スレッド4
列計算スレッド
BLASスレッド
SOTA:「最適化の理論と応用」研究部会
44
dX計算のマルチスレッド化
各列ごとに計算
各列ごとにスレッドで並列計算
今回もBLASのスレッドを無効化
2012/12/15
SOTA:「最適化の理論と応用」研究部会
45
マルチスレッドの効果
SDPA-C
6.2.0(1)
Schur 要素計算
SDPA-C
7.3.8(1)
SDPA-C
7.3.8(2)
SDPA-C
7.3.8(4)
4205.70秒
2094.03秒
1170.37秒
1.78倍
731.98秒
2.86倍
Schur Cholesky
185.36秒
161.03秒
85.24秒
1.88倍
47.77秒
3.37倍
dX の計算
316.00秒
242.51秒
131.03秒
1.85倍
83.54秒
2.90倍
22.22秒
17.93秒
23.88秒
25.86秒
4729.28秒
2515.50秒
1410.52秒
5.31倍高速化
889.15秒
2.82倍
その他
合計時間
1.78倍
（ ）の数字はスレッド数
最大クリーク問題
2012/12/15
SOTA:「最適化の理論と応用」研究部会
46
インデックス
1. 半正定値計画問題の標準形と主双対内点法
2. 行列補完理論と主双対内点法
3. 実装の高速化
1. 行列補完計算の改善
2. Schur complement マルチスレッド化
3. 主探索方向マルチスレッド化
4. 数値実験結果
2012/12/15
SOTA:「最適化の理論と応用」研究部会
47
数値実験環境とテスト問題
CPU Xeon X5365(3.0GHz), Memory 48GB, RedHat Linux
テスト問題1
格子ネットワーク(p,q)で生成したSDP緩和
最大クリーク問題
↑ただし、前処理をする↑
最大カット問題
青が格子の情報で生成される, 赤が変数
テスト問題2
量子化学におけるSpin Glassの基底状態エネルギー計算
DIMACS challenge の torusg-3-15のタイプ
2012/12/15
SOTA:「最適化の理論と応用」研究部会
48
テスト問題の Sparsity Pattern
最大カット問題
クリーク数 1607, 平均8.04, 最大26, 合計 12924
SpinGlass(3D)
クリーク数 591, 平均29.97, 最大773, 合計 17714
2012/12/15
SOTA:「最適化の理論と応用」研究部会
49
最大クリーク問題(p=400,q=10)
クリーク数 438, 平均29.89, 最大59, 合計 13090
SDPA-C
6.2.0
SDPA
7.3.8(4)
SDPA-C
7.3.8(4)
SeDuMi
1.3
Schur
9314.34
262.04
1431.69
-
Cholesky
2601.04
403.28
248.55
-
734.41
13151.10
175.24
-
31.37
12342.98
47.65
-
12681.16
26159.40
1903.13
61861.30
dX
その他
合計
SeDuMi: http://sedumi.ie.lehigh.edu/
SDPA-C が非常に高速に最大クリーク問題を解く
2012/12/15
SOTA:「最適化の理論と応用」研究部会
50
最大カット問題(p=500,q=10)
クリーク数 1607, 平均8.04, 最大26, 合計 12924
SDPA-C
6.2.0
SDPA
7.3.8(4)
SDPA-C
7.3.8(4)
SeDuMi
1.3
Schur
105.25
16.17
317.56
-
Cholesky
502.43
214.43
226.06
-
dX
61.21
3254.05
315.73
-
その他
14.38
3063.60
17.46
-
683.27 6548.26
876.81 95356.20
合計
入力行列が単純なため、マルチスレッドのオーバーヘッドが大きい
SDPA の Schur も非常に高速
2012/12/15
SOTA:「最適化の理論と応用」研究部会
51
SpinGlass(p=15)
クリーク数 591, 平均29.97, 最大773, 合計 17714
SDPA-C SDPA
6.2.0
7.3.8(4)
Schur
SDPA-C SeDuMi
7.3.8(4) 1.3
565.11
7.99
261.40
-
39.47
6.64
11.19
-
dX
496.16
126.90
252.24
-
その他
209.26
194.87
35.570
-
1306.00
336.40
560.40
27898.1
Cholesky
合計
SDPA-C7.3.8は6.2.0よりも高速
ただしSDPAは、さらに高速
2012/12/15
SOTA:「最適化の理論と応用」研究部会
52
SpinGlass(p=25)
クリーク数 3173, 平均29.50, 最大1798, 合計 93609
SDPA-C
6.2.0
Schur
SDPA
7.3.8(4)
SDPA-C
7.3.8(4)
SeDuMi
1.3
35314.96
220.19
11856.50
-
3681.44
520.39
945.49
-
dX
45238.37
11846.59
11594.76
-
その他
23267.71
13436.50
516.45
-
107502.48
26023.67
24913.20
実験なし
Cholesky
合計
このサイズではSDPAよりSDPA-Cのほうが高速
2012/12/15
SOTA:「最適化の理論と応用」研究部会
53
SpinGlass まとめ
計算時間の単位は秒
p
クリーク数
n
平均サイズ
最大サイズ
SDPA7.3.8
SDPA-C7.3.8
10
1000
155
25.69
294
11.85
20.77
15
3375
191
29.97
773
336.40
560.40
18
5832
1118
28.13
913
1522.60
2136.88
20
8000
1737
25.75
1080
3726.03
4552.10
25
15625
3173
29.50
1798
26023.67
24913.20
SDPA-Cのほうが計算時間の上昇が緩やか
大きい問題ほどSDPA-Cが有利
ちなみに、CPU温度は Core i7 3770K で室温約15度の場合
通常 約30度, 最大 約80度, SDPA 約70度, SDPA-C 約60度
2012/12/15
SOTA:「最適化の理論と応用」研究部会
54
まとめ
行列補完理論の数式を改善して高速化
大きい問題ほど効果が高くなる
マルチスレッドでボトルネックを高速化
ただし、性能は疎構造に強く依存
2012/12/15
SOTA:「最適化の理論と応用」研究部会
55