流体アプリの広域分散環境での性能考察

流体アプリの性能考察
遠藤敏夫
背景(1): グリッド計算
グリッド計算の市場規模は2005年には41
億ドル(Grid Technology Partnersの推計)
 実現されているのはサブタスク間の依存
が少ない問題



SETI@home, cell computing etc.
依存関係の多いプログラムはスーパーコ
ンピュータ,クラスタで行われてきた
背景(2): グリッド環境での問題

スパコン・クラスタに比べ,依存関係の多
い問題を解きにくい


計算機が世界中で増減する
通信性能が低い
背景(3): グリッドでの通信性能

一般家庭にブロードバンドが広まればOK?
バンド幅
クラスタ+専用
ネットワーク
クラスタ+普通 グリッド
のネットワーク
>1Gbps
100M~1Gbps
1M~100Mbps
100μs~1ms
10ms~300ms
レイテンシ <10μs

No. バンド幅はクラスタに近づきつつあるが,レ
イテンシはダメ.

地球の裏側(20000km)にはどうやったって66ms以上
(cf. 光速の限界)
本発表の内容

ごく簡単な並列流体シミュレーションプログ
ラムを実装



クラスタでの性能測定


依存関係が多い
Red-black SOR
16CPUで3.5倍(問題サイズ512)
グリッド環境での性能予測

65536CPUで20000倍??(問題サイズ65536)
発表の流れ
流体シミュレーション
 ポアソン方程式の反復解法
 並列アルゴリズムと性能測定
 グリッド環境での性能見積もり
 まとめ

流体シミュレーションとは

流体の挙動を予測
本発表では
 空間は二次元
 各時刻・各位置の速度
を求める
初期条件と境界条件
空間
時間
初期条件

境界条件
与えられた条件(赤)から,未知数(白)を求める
離散化の方法

本発表では,時間幅,空間幅とも一定刻み
として離散化する
有限差分法(今回)
有限要素法
解くべき運動方程式(1)

ナビエ・ストークス方程式
u
1
 (u  )u  p 
u
t
Re
u  0
ただし,u: 速度ベクトル,t: 時間,p: 圧力
Re: レイノルズ数(粘性係数)
  
2
2
   , ,   2  2
x
y
 x y 
解くべき運動方程式(2)

二次元の場合は以下のように式変形できる


    1




t
y x
x y
Re
Ψは流れ関数と呼ばれ,以下が成り立つ


ux 
,uy  
y
x
・・・(#)
プログラムの概要(1)
入力: Ψの初期条件と境界条件
 出力: 各時間,各位置のΨ



Ψが分かれば流速ux, uyはすぐ求まる
データ構造: 二次元配列Ψ, A, B
プログラムの概要(2)
Ψ=初期条件
while (シミュレーション終了) {
for all (x, y) B[x,y]=(#)の右辺;
ΔA=Bを解く
// ポアソン方程式.次で取り上げる
// ここでA=δΨ/δtが成り立つ
for all (x, y) Ψ[x,y]+=A[x,y]*dt;
// Ψ更新(出力)
t += dt;
// 次の時間ステップへ
}
ここでは時間軸方向に前進オイラー法を利用
発表の流れ
流体シミュレーション
 ポアソン方程式の反復解法
 並列アルゴリズムと性能測定
 グリッド環境での性能見積もり
 まとめ

ポアソン方程式とは
 A( x, y )  B( x, y )

の形の方程式
熱拡散,電磁場など応用が広い(らしい)
ポアソン方程式の離散化
( A[ x  1, y ]  A[ x  1, y ]  A[ x, y  1]  A[ x, y  1]  4 A[ x, y ]) / h 2
 B[ x, y ]
・・・(*)
が全点で成り立つ (hはメッシュ幅)
これは疎行列連立方程式
であり,反復法で解くのが普通
h
反復法とは
for all (x,y) 適当にA[x,y]を決める
while (true) {
for all (x,y) {
A’[x,y] = 「A[x,y]よりも良く(*)を満たしそうな値」
}
if (全点において |A’[x,y]-A[x,y]|<ε) break;
// 収束した!!
for all (x,y) A[x,y] = A’[x,y]
// 次回のため
}
※ 収束までに必要な回数は方式・問題によって異なる
反復法(1):ヤコビ法


どうやって良いA’[x,y]を見積もるか?
(*)を変形し,A[x,y]をA’[x,y]で置き換える
A'[ x, y ] 
( A[ x  1, y ]  A[ x  1, y ]  A[ x, y  1]  A[ x, y  1]  h 2 B[ x, y ]) / 4
欠点:収束が遅い.二次元メッシュの場合はO(n2)
回の反復
 周りの4点に古い値を使っているため
反復法(2): Red-Black法


メッシュ上の点を交互に塗り分
ける
各反復において,



赤の点をヤコビ法と同様に更新
黒の点を,新しい赤の点を用いて
更新
ヤコビ法より収束が速いものの,
依然O(n2)回
反復法(3): SOR法
従来の方式
A[x,y]



SOR
A[x,y]
従来の方式通りに計算し,更新時の変更を定数
ω倍とする(1<ω<2)
最適なωのとき,反復回数はO(n)に
以下の評価ではRed-black SORを並列化
反復法(4):マルチグリッド法



本来の
メッシュ

粒度の異なるメッシュを
用意
粗いメッシュで大まかに
計算し,結果をより細か
いメッシュへ反映
収束はとても速い
SPLASH2ベンチマーク
のocean
発表の流れ
流体シミュレーション
 ポアソン方程式の反復解法
 並列アルゴリズムと性能測定
 グリッド環境での性能見積もり
 まとめ

Red-black SORの並列化

メッシュをプロセッサたちで分担

プロセッサ境界を越える依存関係があるの
で,メッセージ通信が必要
並列Red-black SORの通信パ
ターン
各反復において,
 赤更新の前
 黒更新の前
に通信が必要
クラスタにおける性能

Exec. time (s)
100
Sun Bladeクラスタ

80

60
40
UltraSPARC 750MHz
100Mbps ether

Phoenixライブラリ利用

16プロセッサで2.8倍の
台数効果.遅い!
通信に80%の時間を使う
プロセッサあり
20
0
0
5
10
15
# of processors
1タイムステップの計算時間
メッシュサイズ512×512
反復回数2072となった
20

通信をさぼって速くする工夫(1)
複数列(f列)を一度に送ってしまう
→ レイテンシの影響を1/fにできる

通信
計算
時間
通信をさぼって速くする工夫(2)
欠点:
 だぶった領域の計算が必要で,計算コストが増える
 バンド幅利用は不変
計算すべき領域
Exec. time (s)
工夫後の性能
120

100

80
先ほどより向上
f=8のとき,16プロセッ
サで3.5倍の性能
60
40

20
0
0
5
10
15
# of processors
FETCH-1
FETCH-4
FETCH-2
FETCH-8
20
まだまだ解析・改善の
必要
発表の流れ
流体シミュレーション
 ポアソン方程式の反復解法
 並列アルゴリズムと性能測定
 グリッド環境での性能見積もり
 まとめ

果たしてグリッド環境で性能は
出るのか?

広域で,多数のコンピュータを使った実験
は必要な予算・手間が多大

苦労の後に「やっぱり性能出ませんでした」と
いうのは悲しい
→性能の見積もりが必要
メッシュサイズn×nのとき,
pプロセッサでの性能は?
仮定
見積もりを楽にするために以下の仮定を置く
 計算速度は一定値S (1データ,1回の更新あたり)



レイテンシは一定値L



キャッシュなどの影響なしとする
Sunクラスタでの性能を基に,S=0.2usとする
L=100msとする
バンド幅は無限大
各プロセッサの担当領域は正方形
一反復あたりの実行時間の見
積もり



実行時間=通信時間+計算時間
通信時間=2L/f
計算時間=S(b+2f-2)2
ただし,
b: 各プロセッサが担当する領域の一辺.b=n/√p
f: 同時に通信する列数
性能の見積もり結果(1)
通信数削減なし(f=1)の場合の台数効果

右のグラフは対数軸
60000
100000
50000
10000
40000
1000
speed-up
speed-up

30000
20000
10000
100
10
1
0
0.1
0
n=1024
n=65536
20000
40000
60000
# of processors
n=4096
n=262144
80000
n=16384
1
n=1024
n=65536
10
100 1000
# of processors
n=4096
n=262144
10000 100000
n=16384
性能の見積もり結果(2)

メッシュサイズが大きいほどスケーラビリティ良


通信数に比べ計算量の比率が高いため
十分なスケーラビリティを得るにはn≧256K
64kプロセッサで
 n=256k → 28000倍
 n=64k → 3000倍
性能の見積もり結果(3)

通信数削減ありの場合の台数効果
最も実行時間(2L/f+S(b+2f-2)2)を小さくするようなfを
用いるとする
60000
100000
50000
10000
40000
speed-up
speed-up

30000
20000
1000
100
10000
10
0
1
0
n=1024
n=65536
20000
40000
60000
# of processors
n=4096
n=262144
80000
n=16384
1
n=1024
n=65536
10
100 1000
# of processors
n=4096
n=262144
10000 100000
n=16384
性能の見積もり結果(4)
レイテンシの影響を小さくすると,何が変わった?
 より小さい問題サイズでも,性能を絞りだせる!
64kプロセッサで
 n=64k → 22000倍
発表の流れ
流体シミュレーション
 ポアソン方程式の反復解法
 並列アルゴリズムと性能測定
 グリッド環境での性能見積もり
 まとめ

まとめ

Red-black SORを用いた流体アプリについて,



クラスタ上での性能測定
グリッド環境での性能予測
CPUコストを増やしてでも,通信レイテンシの影
響を少なくする方向は正しいと思われる
課題

性能モデルの改善



バンド幅を考えないのはあまりに非現実的
遅延のばらつきに適応可能なアルゴリズムの考
案
遅延の影響を吸収できるようなタスクスケジュー
ラ?
参考文献


J. P. Singh et al., “SPLASH: Stanford Parallel
Applications for Shared-Memory”, In
Computer Architecture News, vol. 20, no. 1.
William Press et al.,「NUMERICAL RECIPES in
C,日本語版」,技術評論社,1993