並列計算法 （Parallel Computing） スーパーコンピュータの性能 並列

計算流体力学
（第13回資料）
2014年1月14日
スーパーコンピュータの性能
33862Tflops(2013)：Tianhe-2
並列計算法
（Parallel Computing）
Flops
1Pflops 10510Tflops(2011)：京
36Tflops(2002)：地球シミュレータ
1Tflops
First Parallel Computer
1972
話のポイント：
・並列計算機の仕組みと並列化の方法
Massive parallel
multi-processor
vector
並列計算の目的
・より大規模な問題の計算が可能となる．
(n倍のメモリー，n:プロセッサ数)
・計算時間を大幅に短縮できる．（計算時間1/n）
scaler
micro processor
vacuum tube
transistor
http://www.top500.org
並列コンピュータの必要性
並列計算の歴史
科学的に複雑な問題の解明
１９８０年代末に提案されたグランド・チャレンジ問題
・計算流体力学と乱流
・大域的な気象と環境のモデル化
・物質の設計と超流動
・医学，人間の臓器と骨格のモデル化
・宇宙論と宇宙物理学
：
コンピュータの性能向上
回路素子の進歩
真空管→トランジスタ→集積回路→超大規模集積回路
Richardson（英国）
１９１１年：非線形偏微分方程式の数値解法を発表
１９２２年：「数値的方法による天気予測」
・大気を５層に分け，水平方向にはヨーロッパ大陸を
ブロックに分割した
・６時間分の計算を手回し計算機で６週間
シングルプロセッサの限界
真空中の光の速度（約 3 108 m/sec)
Richardsonの夢
北半球全体を約２，０００個のブロックに分け，３２人が
一つのブロックの計算を担当すれば，６時間先の予報を
３時間で出来ると見積もった．
（６４，０００人が巨大な劇場に集合し，指揮者のもとで
一斉に計算を行えば実現できる）
Richardsonの夢
Richardson,L.F.:Weather Prediction by Numerical Process,
Cambridge University Press, London, 1922
並列コンピュータの出現
イリノイ大学Daniel Slotnic：２種類の並列コンピュータを設計
１９７２年：世界最初の並列計算機ILLIAC Ⅳ(Burroughs社）
完成（６４台の処理ユニット，１台の中央処理ユニット：SIMD)
日本での並列コンピュータ構築研究
１９７７年：PACS/PAXプロジェクトが発足
星野力（京都大学→筑波大学）
・PACS-9(１９７８年，０．０１Mflops）
・PAX-32(１９８０年，０．５Mflops）
１９８０年PACS/PAXプロジェクト筑波大学
に移動
・PAX-128（１９８３年，４Mflops)
・ PAX-64J（１９８６年，３．２Mflops)
・QCDPAX(１９８９年，１４Gflops)
・CP-PACS（１９９６年，３００Gflops）
（１９９７年，６００ Gflops：２０４８CPU）
地球シミュレータ（初代）
地球シミュレータ（現在）
640 nodes
（5120 CPU）
8GFLOPS@5120
＝40TFLOPS
35.86Tflops
(実効ピーク性能)
http://www.jamstec.go.jp/esc/
２００９年３月導入の新たなシステム（NEC SX-9/E）
2002年から2年半世界最高速
次世代スパコン 京
世界472位（2013年11月現在）
160 nodes（1280 CPU）
102.4GFLOPS@1280＝131TFLOPS
122.4Tflops(実効ピーク性能)
TSUBAME 2.5(東工大)世界11位
2843TFLOS（実効ピーク）
共有メモリー型と分散メモリー型
共有メモリー型（shared memory type）：すべてのプロセッサー
がメモリーを共同で使用するタイプ
分散メモリー型（distributed memory type）：すべてのプロセッ
サーが個別にメモリーを持っている
世界4位（2013年11月現在）
http://www.nsc.riken.jp/index_j.html
705024 cores（88128 CPU）
128GFLOPS@88128＝11280TFLOPS
10510Tflops(実効ピーク性能) 実効効率93％
並列計算機とプログラミングの方法
並列計算機とプログラミングの方法
共有メモリシステム
共有分散メモリシステム（shared-distributed memory type）:
大規模並列システムの場合
OpenMP
MPI
MPI
分散メモリシステム
MPI
OpenMP
MPI
PC cluster
OpenMP
MPI
OpenMP
MPI
GPUを用いた汎用計算
• GPGPUとは
– General Purpose Graphics Processing Unit
– GPUを一般的な用途で利用することの総称
PC(Gateway)
• GPUの特徴
HUB
– 演算能力が非常に高い
• 例：GeForce GTX Titan：4.5TFlops
OS:Linuxが主流
PCs
自作も容易
– デスクトップPC(Win, Mac & Linux)に搭載可能．
（もちろんノートPCにも）
– 高FLOPSの割に値段が安い
• GeForce GTX780(4.0Tflops)：50,000円程度
• 数値計算専用のモデルもある
（Tesla シリーズ，搭載メモリ多い) Tesla K40(スパコン向け)
GPU計算への応用(CUDA FORTRANを利用)
SM:ストリーミング
マルチプロセッサ
448/32=14 個
SP:ストリーミングプロ
セッサ
SM1つあたり32個
・ハードウエア
従来のベクトルスーパーコンピュータの性能をはるかに
上回る性能の実現（マイクロプロセッサの性能向上，高速
ネットワーク）
超高性能並列コンピュータからPCクラスターまで
シェアードメモリ
(超高速)
レジスタ・メモリ(超高
速)
＊NVIDIA製GPU
のアーキテクチャ
144GB/sec
並列計算機の性能向上と普及
・ソフトウエア
MPI，OpenMPなどの並列プログラミング支援ツールの整備
ビデオメモリ（グローバルメモリ）
プロセッサ
SP数
クロック周波数
ピーク演算性能
メモリ 17
標準メモリ設定
メモリバンド幅
Tesla C2075(GPU)
448個
1147MHz
1027.7GFLOPS(単精度)
515GFLOPS(倍精度)
6GB GDDR5
144GB/s
Intel Core i7 990Xの場合
← 6個
← 3460MHz
← 83.2GFLOPS
誰もが並列計算できる環境が整った（パラダイムの変化）
欧米のほとんどの土木・機械系の大学院では並列計算
に関する講義科目が用意されている
← 25.6GB/s
並列プログラミング支援ツール
１）OpenMPを用いる（共有メモリ型）
プログラム構造を変更することなく指示文を挿入するだけで
並列化可能
２）メッセージパッシングライブラリーを用いる
（分散メモリ型，共有メモリ型）
MPI(Message Passing Interface),
ライブラリー関数群（逐次言語（Fortran,Cなど）の
ライブラリーとして使用)
本格的な並列プログラミングが可能
３）自動並列化コンパイラーを用いる（共有メモリ）
プログラム（逐次言語）をコンパイラーが自動並列化
OpenMPとは
・共有メモリ型並列計算機に対する並列化手法
・コンパイラに対して並列化の指示を行なうために，プログラム
中に追加する指示行の記述方法を規定
（ FortranやC/C++の並列化に対する標準規格）
特徴：
１）並列化が容易
逐次プログラムに指示行を挿入
２）移植性が高い
OpenMPに準拠しているコンパイラ
３）逐次プログラムとの共存が容易
並列プログラムとしても逐次プログラムとしても利用できる
OpenMPによる並列計算
OpenMPプログラムの並列計算の概要
マスタースレッド
OpenMPの基本的な規則
!$OMP
○○○
!$
○
!
○×△
××× ：OpenMP 指示文
= ×＊△
：OpenMP のみで実行される演算．
program welcome_to_parallel_world_omp
integer :: mythread
mythread = 0
write(6,*) mythread, ‘Hello Parallel World’
write(6,*) ‘Hello Parallel World’
：コメント文．
!$OMP parallel private( mythread )
program hello_world_omp
integer :: mythread
mythread = 0
!$OMP parallel private( mythread )
!$ mythread = omp_get_thread_num()
wirte(6,*) mythread,‘Hello World’
!$OMP end parallel
end program hello_world_omp
Parallel構文
(並列リージョンの立上げ)
!$ mythread = omp_get_thread_num()
write(6,*) mythread,‘Parallel Computing’
!$OMP end parallel
環境変数の獲得
（自分が何番のスレッドなのか）
Parallel構文
(並列リージョンの終了)
メッセージ通信によるプログラミング
プロセッサ間やネットワークで結ばれたコンピュータ間の制御，
分散主記憶に分散された情報の交換を，メッセージ通信により
行なう．
分散メモリー型並列計算機では，メッセージ通信によるプロ
グラミングが最もハードウエア構成に即したプログラミングで
あり高性能な計算が行える
（共有メモリー型並列計算機でも有効）．
MPI(Message Passing Interface)が標準
（MPI-1:基本仕様，MPI-2:拡張機能仕様)
無償でダウンロードできる
スレーブスレッド
write(6,*) mythread,‘Good-by’
end program welcome_to_parallel_world_omp
write(6,*) ‘Good-by’
MPI（Message Passing Interface）とは
・ノード間の情報交換をメッセージ通信で行なうライブラリ
（プログラム言語ではない）
・FortranやCなどの従来の逐次プログラム言語中において，
ライブラリとして引用する
（各ノードにおいて不足する情報の交換を，MPIサブルーチンを
適切にcallすることにより行なう）
・様々なプラットホームで利用できるフリーウェア
MPIの実装
MPICH：http://www-unix.mcs.anl.gov/mpi/mpich/
LAM:http://www.lam-mpi.org/
Windows環境で動作するMPIもある
並列計算プログラミング
MPIのサブルーチン
１）環境管理サブルーチン
MPI_INIT, MPI_COMM_SIZE, MPI_COMM_RANK
MPI_FINALIZE
例題
２）グループ通信サブルーチン
MPI_REDUCE, MPI_ALLREDUCE
配列A
３）1対1通信サブルーチン
MPI_ISEND, MPI_IRECV, MPI_WAITALL, MPI_WAIT
４）その他
MPI_WTIME
：
配列の和を求める
【 配列A(1)からA(100)の総和を求める 】
A(1) A(2) A(3)
……………………………………………… A(100)
逐次プログラム
do i=1,num
a(i)=i
end do
複数のプロセッサで
並列に処理し、
実行時間の短縮を図る
s=0
do i=1,num
s=s+a(i)
end do
並列計算プログラミング
implicit double precision (a-h,o-z)
parameter ( num = 100 )
integer a(num)
integer s
プログラム名：single_sum.f
並列計算プログラミング例
配列をプロセッサに分割する
配列A(100)の和の計算を４プロセッサ(PE)に分担させた。
A(1)+
+ A(100)
c
do i=1,num
a(i)=i
end do
s=0
do i=1,num
s=s+a(i)
end do
c
write(6,100) s
100 format(3x,'sum = ', i10)
c
stop
end
A(1)………A(25)
PE 0 が担当
A(26)……A(50)
PE 1 が担当
A(51)……A(75)
PE 2 が担当
部分和 + 総和
s = A(1)+…+ A(25)
部分和
s = A(26)+…+ A(50)
部分和
s = A(51)+…+ A(75)
PE:0
PE:1
PE:2
PE:3
CPU
CPU
CPU
CPU
メモリー
A(100), S
メモリー
A(100), S
メモリー
A(100),S,TS
メモリー
A(100), S
A(76)……A(100)
PE 3 が担当
部分和
s = A(76)+…+ A(100)
総和 TS= S(PE0) + S(PE1) + S(PE2) + S(PE3)
並列計算プログラミング例
MPIプログラムの基本的構造
変数の定義,宣言
include 'mpif.h'
c
→MPIプログラムのコンパイルに必要
call MPI_INIT( ierr )
c
→MPI環境の初期化
call MPI_COMM_RANK ( MPI_COMM_WORLD, irank, ierr )
c
→各プロセスのランクirankを取得
call MPI_COMM_SIZE ( MPI_COMM_WORLD, isize, ierr )
c
→プロセス数isizeを取得
:
並列処理プログラム
:
call MPI_FINALIZE( ierr )
c
→MPIの終了
stop
end
並列計算プログラミング例
do i = ist, iet
S = S + A(i)
end do
c
integer A( num )
integer S
integer TS
c
include 'mpif.h' →MPIプログラムのコンパイルに必要
call MPI_INIT( ierr )
→MPI環境の初期化
call MPI_COMM_RANK ( MPI_COMM_WORLD, irank, ierr )
→自分のランクirankを取得
call MPI_COMM_SIZE ( MPI_COMM_WORLD, isize, ierr ) →プロセス数isizeを取得
c
do i = 1, num
A(i) = i
end do
c
S=0
TS = 0
ipe = num / icp
ist = 1 + irank * ipe
iet = ( irank + 1 ) * ipe
MPIサブルーチン引用例
call MPI_REDUCE( S, TS, 1, MPI_INTEGER, MPI_SUM, 0,
&
MPI_COMM_WORLD, ierr )
→各プロセスの部分和Sをプロセス０（ランク０）に集め総和を求めTSに格納
c
call MPI_REDUCE( S, TS, 1, MPI_INTEGER, MPI_SUM, 0,
&
MPI_COMM_WORLD, ierr )
→各プロセスの部分和Sをプロセッサ０に集め総和を求めTSに格納
c
write(6,100) irank, S, TS
c
100 format(3x,'irank = ',i3,3x,2i10)
c
call MPI_FINALIZE( ierr )
→MPIの終了
c
stop
end
implicit double precision (a-h,o-z)
parameter ( icp = 4 )
parameter ( num = 100 )
計算結果
irank =
irank =
irank =
irank =
0
1
2
3
325 5050
950 0
1575 0
2200 0
call MPI_REDUCE(sendbuf, recvbuf, count, datatype, op, root,
comm, ierror)
・sendbuf:送信バッファの先頭アドレスを指定する.
・recvbuf：受信バッファの先頭アドレスを指定する.
・count：整数.送（受）信メッセージの要素数を指定する.
・datatype：整数.送（受）信メッセージのデータタイプを指定する.
・op：整数.演算の種類を指定する.
・root：整数.宛先プロセスのcomm内でのランクを指定する.
・comm：整数.送受信に参加するすべてのプロセスを含むグループの
コミュニケータを指定する.
・ierror：整数.完了コードが戻る.
連立一次方程式の並列解法
陽解法
並列化の容易さ：容易
並列化効率：良好
陰解法
直接法
並列化の容易さ：やや難
並列化効率：やや不良（超並列には不向き）
並列プログラミング参考書
・樫山，西村，牛島：並列計算法入門，日本計算工学会編
（計算力学レクチャーシリーズ３），丸善（2003）
１章
２章
３章
４章
５章
付録
計算力学における並列計算概論
有限要素法における並列計算法
有限差分法における並列計算法
境界要素法における並列計算法
PCクラスターの構築法とグリッドコンピューティング
並列プログラム
Melosh et al., Comp. Struct., Vol.20, pp.99-105, 1985
ガウスの消去法が最も効率がよい．
反復法
並列化の容易さ：容易
並列化効果：良好
Hughesらによって提案されたelement-by-element法に基づく
前処理付き共役勾配法…省メモリーでかつ高効率
Barragy and Carey, Int.J.Numer.Meth.Engng., Vol.26, pp.2367-2382, 1988
Tezduyar et al.1992-present（流れの並列計算に関する論文多数）
並列計算効率の評価指標
・牛島：OpenMPによる並列プログラミングと数値計算法，丸善
（2006）
1章
2章
3章
4章
OpenMPと並列計算
OpenMPによるプログラミング
OpenMPによる並列数値計算
付録
並列計算における計算時間
全計算時間＝計算時間＋通信時間＋同期待ち
スピードアップ =
並列化効率 =
１ノード使用時の実行時間
N ノード使用時の実行時間
スピードアップ
N
粒度(granularity)：並列処理の単位
粗粒度：通信負荷に比べて計算負荷が卓越
全計算時間は遅いが並列化効率は良い
細粒度：通信負荷が増大
全計算時間は速いが並列化効率は悪い
同期待ち 時間
実行時間＝計算時間＋通信時間＋同期待ち時間
計 算時間
: :
並列計算における実行時間
ノード １
通信時間
ノードｎ
実行時間
1ス テッ プあたりの 実行時間
粒度が小さい（細粒度）と並列化効率が悪くなる理由
高効率な並列計算を実現するには
１）並列処理可能な部分の拡大
可能な限りの処理の並列化
実行時間
計算時間（並列化不可能部分）
ｎ/ 1
計算時間（並列化可能部分）
通信時間
ｎ台
1
台
並列化可能
並列化不可能
計 算 時間
２）各プロセッサの計算負荷の均等化
各プロセッサーに割り当てる要素数の均等化
３）プロセッサ間の通信量の最小化
1
2
4
8
1 6 ノード数
通信を行う節点の最小化
自動領域分割法が有効
アムダールの法則
P：並列化部分率, N：使用するプロセッサ数
1
スピードアップ 
(1  P )  P / N
領域分割法
MPIに基づく並列計算では領域分割に基づく並列計算が一般的
解析領域
有限要素方程式
小領域の境界
自動領域分割法
(1) 小領域内での要素数の均等化
(2) 小領域間の境界上の節点数の最小化
効率の良い並列計算
伊勢湾台風による被害 (1959)
領域分割法
（a）節点ベース
（b）要素ベース
(a) 節点ベース
9
10
11
(b) 要素ベース
12
9
10
11
12
5
6
7
8
1
2
3
・1
5
6
7
1
2
3
8
・2
4
4
要素ベースの方法
Bisection法：計算領域を次々に数学的手法を用いて2分していく
方法．
Greedy法：隣接要素を記憶しておき近隣の要素を必要数だけ
探索する方法．数学理論を用いない簡便な方法
有限要素メッシュ
（elements:206,977，nodes:106,577）
領域分割図
（512 sub-domain）
領域分割に基づく並列計算例
領域分割に基づく並列計算例
計算結果と観測結果との比較（名古屋港）
速度向上比（Speed-up ratio）
速度向上比 =
１プロセッサ使用時のCPU時間
N プロセッサ使用時のCPU時間
並列化効率（Efficiency）
並列化効率=
速度向上比
N
領域分割図
大気環境流れ
通信ライブラリ ： MPI
自動領域分割システム ：Metis
例：64分割図
まとめ
並列計算機環境の進歩・普及と並列化手法ついて紹介した．
・マルチコアCPU登載のパソコンの普及により、誰もが並列計算
可能な環境にある
・共有メモリー型コンピュータの場合には、OpenMPが、分散メ
モリー型コンピュータの場合にはMPIが有効である
・大規模計算には領域分割に基づく並列計算（MPI）が並列化効率
の点から一般的である
・連立一次方程式の解法は反復法が並列化効率の点で優れる
今後の期待
・ハイブリッド計算（Open MP+MPI）
・GPUコンピューティングのさらなる進歩・普及
商用ソフトウエアの紹介
世界でもっとも使用されているソフトウエア
1位 FLUENT
有限体積法に基づく使用実績世界No.1の汎用熱流体解析ソフト
ウェア．
http://ansys.jp/products/fluent/
2位 STAR-CD（日本では最も使用されている）
有限体積法に基づく汎用熱流体解析ソフトウェア．
http://www.cdaj.co.jp/
http://acusolve.jsol.co.jp/index.html
自動車系に強い
AcuSolve：有限要素法（Galerkin/Least-Squares：GLS法）による
汎用熱流体解析ソフトウェア
http://acusolve.jsol.co.jp/index.html
AcuSolve
AcuSolveは有限要素法による汎用熱流体解析システム
・有限要素法（Galerkin/Least-Squares：GLS）法による汎用性と高精度な解
・テトラメッシュの使用を推奨
・メッシュの歪みに頑強で、実用形状のメッシュ作成工数を大幅に低減できる
・反復収束が速やかで解析時間を予測できるので、製品開発設計工程での
計画的な適用が可能
・３次元CADとの連携やユーザーインターフェースのカスタマイズによって
CAEオートメーションが可能
http://acusolve.jsol.co.jp/index.html
STAR-CD
STAR-CD英国インペリアル大学のゴスマン教授を中心として開発された有限体積法
非構造メッシュによる汎用熱流体解析プログラム．http://www.cdaj.co.jp/
CD-adapcoにより開発
STAR-CDがカバーする解析対象は単純な熱流れ解析から，固体熱伝導や輻射／
太陽輻射を含んだ熱問題，様々なタイプの混相流問題，化学反応／
燃焼問題，回転機器問題，流体騒音問題，移動境界問題など多岐に渡り，
解析メッシュ作成機能としてもサーフェスメッシャー，サーフェスラッパー，
ボリュームメッシャーなどを標準機能として提供するほか，CADとの親和性を高めた，
CADの中でSTAR-CDを利用するアドオン環境（CATIA v5, Pro/E, UG-NX, SolidWorks）
も用意．
非構造／完全不連続メッシュ（メッシュのつながりを意識せずにモデル化が可能）
による複雑形状への優れた柔軟性とメッシュ作成コストの削減（ヘキサ，テトラ，
プリズム，トリムセルを含む多面体セルに対応）
自動車，航空宇宙，重工重電，家電，化学，建築，官公庁など
（国内1200ライセンス，全世界で3500ライセンスを超える導入実績）
FLUENT
・FLUENTは、有限体積法に基づく使用実績世界No.1のCFD(数値流体力学)ソフトウェア。
http://ansys.jp/products/fluent/
・FLUENTには豊富な物理モデルが搭載されており、化学反応、燃焼、混相流、
相変化などが取扱える。また、高度なカスタマイズも可能。
・プリプロセッサには、GAMBITを標準装備。
FLUENTでは連成解法および分離解法のアルゴリズムを用いて、
非構造格子で離散化した質量、運動量、熱、化学種に関する
支配方程式を解く。
・2次元平面、2次元軸対称、旋回を伴う2次元軸対称、および3次元の流れ・定常/非定常流
・あらゆる速度域に対応 ・非粘性流、層流、乱流 ・ニュートン流体、非ニュートン流体
・強制/自然/混合対流を含む熱伝達、固体/流体連成熱伝達、ふく射・化学種の混合/反応
・自由表面流/混相流モデル・分散相のラグランジュ式追跡計算
・融解/凝固を対象とした相変化モデル ・非等方浸透性、慣性抵抗、固体熱伝導、
空隙率を考慮した速度、多孔質面圧力ジャンプ条件を含む多孔質モデル
・ファン、ラジエータ、熱交換器を対象とした集中定数モデル
・移動物体周りの流れをモデル化するダイナミックメッシュ機能 ・慣性または非慣性座標系
・複数基準座標系およびスライディングメッシュ
・動静翼相互作用をモデル化する接続境界面モデル
・流体騒音予測のための音響モデル ・質量･運動量･熱･化学種の体積ソース項
・物性値のデータベース