反応・相変化を伴う多分散系混相流シミュレーションコードの最適化
佐々木季広 2)，磯部洋子 2)，弘貴 3)，渡部広吾輝 3)，志村啓 3)，奥野亜紀子 3)，江川隆輔 4), 5)，滝沢寛之 4), 5)，小林笠原松尾大輔 1)，加藤航平 3)，広明 5), 4)1) 東北大学 情報部情報基盤課
2) 日本電気株式会社
3) 慶應義塾大学
4) 東北大学 サイバーサイエンスセンター
5) 東北大学 大学院情報科学研究科
d-sasaki@cc.tohoku.ac.jp

Performance Optimization of Large-scale Simulation of Polydisperse
Multiphase Flow with Reaction and Retastability
Daisuke Sasaki1), Toshihiro Kato2), Yoko Isobe2), Hirotaka Kasahara3), Hiroaki Watanabe3),
Kei Shimura3), Kohei Okuno3), Akiko Matsuo3), Ryusuke Egawa4), 5),
Hiroyuki Takizawa4), 5), Hiroaki Kobayashi 5), 4)
1) Information Infrastructure Division of Information Department, Tohoku Univ.
2) NEC Corporation.
3) Keio University.
4) Cyberscience Center, Tohoku Univ.
5) Graduate School of Information Sciences, Tohoku Univ.
概要
東北大学サイバーサイエンスセンターでは，大規模科学計算システムを最大限活用するために計算科学者と本センターの計算機科学者が連携しながら，プログラムの高速化技法の研究・開発に取り組んでいる．本稿では，慶應義塾大学の研究グループが研究・開発している反応・相変化を伴う多分散系混相流シミュレーションコードのベクトル型スーパーコンピュータ SX-ACE における高速化について紹介する.スーパーコンピュータを利用することで，ワークステーションでの実行よりも圧倒的に大規模なモデルを高速に解析することが可能となり，粉塵爆発現象に対する安全防護技術の発展に貢献することが期待される.

図 1 スーパーコンピュータ SX-ACE

1

はじめにリバンド幅をそのコアで占有して 利用でき，4Bytes/FLOP での実行が可能となっている．SX-ACE

東北大学サイバーサイエンスセンター(以下，本は高ベクトル演算性能と，高データ供給性能によセンター)は，全国共同利用機関として大規模科学計算システムの運用と，本システムを最大限に活り，高い実行性能を達成することができる．
本センターの SX-ACE は 2,560 ノードで構成さ用可能なコードの高速化技法や新しいシミュレーれている．各ノード間は最大 4 GB/sec×2(双方向)ション技術の研究・開発に取り組んでいる．本センターでは 1997 年からさまざまな分野における実アプリケーションの高速化支援を計算科学分野で接続された 2 段ファットツリーネットワークで構成され，システム全体の性能は 706.6 TFLOPS となっている．利用者には，1 ジョブあたり最大1,024 ノードの大規模な実行環境を提供している．の利用者との共同研究を通じて行っている[1]．また，センター独自の共同研究に加え，全国の情報基盤センター等と連携して，学際大規模情報基盤表 1 SX-ACE の諸元CPUあたりコア数
理論最大演算性能
最大ベクトル演算性能
ADB
メモリバンド幅(コア-ADB間）
メモリバンド幅(CPU-メモリ間）

4個
276GFLOPS
256GGLOPS
1,024KB×4
256GB/s×4
256GB/s

ノードあたり

CPU数
メモリ容量

1個
64GB

システムあたり

ノード数
CPU数
メモリ容量
理論最大演算性能
最大ベクトル演算性能

2,560ノード
2,560個160TB
706.6TFLOPS
655.4TFLOPS

共同利用・共同研究拠点(JHPCN)や革新的ハイパフォーマンス・コンピューティング・インフラ(HPCI)を構成し，多様なニーズに応える計算環境の提供も行っている．以上のように本センターでは，利用者である計算科学者と本センターの計算機科学の専門家が密接に連携しながら，科学・工学の恒常的な進歩を支える共同研究，利用者コードの高速化支援活動を推進している．
本稿では，慶應義塾大学の研究グループが研究・開発している反応相変化を伴う多分散系混相流シミュレーションコードの高速化について紹介する．

3 反応・相変化を伴う多分散系混相流シ
ミュレーションコード
本コードでシミュレーションを行っている粉2

スーパーコンピュータ SX-ACE

塵爆発とは，可燃性の固体粒子が一定条件下で空気中に浮遊している状態で，何らかの現象により本センターに導入されているスーパーコンピュ発火し爆発する現象のことである．工場等のもの

ータ SX-ACE を図１に，諸元を表 1 に示す．SX-ACE

つくり現場において，大規模な爆発を伴う火災がの CPU はこれまでのベクトルプロセッサと同様に，高ベクトル演算性能と高メモリバンド幅を継承し発生する危険がある.この爆発を伴う火災は，工場内で発生した大量の粉塵が浮遊・堆積し，静電気などにより発生した火花に引火することが原因でている．また，CPU は 4 コアで構成され最大ベク発生する.また，この火災で生じる衝撃波により堆トル演算性能は 256 GFLOPS である．各コアには容積していた粉塵が舞い上がりさらに大規模な火災量が 1,024 KB の ADB(Assignable Data Buffer)がに発展する．

搭 載 さ れ ，ADB とコ ア 間の メ モ リバ ン ド幅 は256GB/sec を有する．データが ADB 経由でアクセスされる場合は， データ転送速度と計算速度の比

慶應義塾大学の研究グループでは，粉塵爆発現象について，粉塵からの可燃性ガスの揮発及び炭素個体の燃焼や気相燃焼，揮発ガス・既燃ガス・空気の混合を考慮しシミュレーションを行っていである B/F 値が 4 Bytes/FLOP となり，メモリ負荷る．また，爆発に伴って生じる衝撃波，爆風によの高いアプリケーションでも高い実行効率で実行る粉塵の巻き上げについて解明するため，シミュ可能となっている．また，単一コアで演算を実行レーションコードの開発に取り組んでいる．

する場合には，CPU-メモリ間の 256 GB/sec のメモ

４ コードの高速化
4.1 SX-ACE における初期性能
慶應義塾大学の研究グループが開発したプログ
ラムの初期性能解析結果を表 2 に示す．性能解析には SX-ACE の簡易性能解析ツール FTRACE を用いた．また，測定データのグリッドサイズは 1,901×161 である．
初期性能解析の結果に基づき，全体の実行時間に対して約 90%を占めているサブルーチン A～C に

図 2 サブルーチン A,B のオリジナルコードついて高速化を行うことを検討する．また，将来の大規模並列実行を見据え，コードの並列化も検両サブルーチンの最内ループは 3 重ループの最内ループとなっているため，最内ループの外側の討する．

ループでベクトル化することで性能の改善が期待
表2

初期性能解析

ベクトル ベクトル化率
長subroutine A
8039.993
91.75
21.9
subroutine B
2314.420
82.60
18.8
subroutine C1164.311
13.23
12.3
other subroutines 1850.370
total13397.172
88.92
26.2

サブルーチン名 実行時間できる．今回は図 3 のように unroll 指示行を挿入し，最内ループを展開，最内ループの外側のループでベクトル化することで性能改善を試みる．

以下では，本コードに適用した高速化について述べる．
4.2 サブルーチン A，B
表 2 からサブルーチン A，B はベクトル長が非常に短いことが分かる．ベクトルプロセッサは複数のデータに対して一括で演算を行うことができ，

図 3 サブルーチン A,B の高速化後のコード最内ループの反復(イテレーション)回数であるベ
ク トル 長が 長い ほど 演 算の 効率 が高く なる ．

高速化後の性能を表 3 に示す．指示行を挿入す

SX-ACE ではベクトル長を 256 の時に演算効率が最ることでベクトル長がサブルーチン A で 21.9 から大になる．

225.9 に，サブルーチン B で 18.8 から 250.9 まで

サブルーチン A，B の内部を詳細に解析したとこ
ろ，両サブルーチンは同様のループ構造となって大きく改善され，サブルーチン A で約 44 倍，サブルーチン B で約 17 倍の性能向上を得た．

おり，非効率なベクトル化が行われている箇所があった(図 2)．
このループの反復回数は 10 であり，表 3 高速化の効果ベクトル化されているものの，SX-ACE の理想のベ実行時間 ベクトル ベクトル 加速化率
長率
[sec]
before 8039.993
91.75
21.9
subroutine A
44.17
after
182.020
99.18
225.9
before 2314.420
82.60
18.8
subroutine B
16.96
after
136.431
99.97
250.9

クトル長より非常に短く，ベクトル化の恩恵を受けていない．

4.3 サブルーチン C

表 2 からサブルーチン C はベクトル化率が低く，ベクトル長が短いことが分かる．サブルーチン C の内部を詳細に解析したところ，ループ内に依存関係があり，最内ループがベクトル化されていないことが分かった(図 4）．また，最内ループは外側のループ毎にループ長が異なり，かつそのループ長は短い．外側ループのループ長が長いため，外側ループを最内ループに移動することで性能改善が期待できる．

図 5 サブルーチン C のベクトル化後のコード

図 4 サブルーチン C のオリジナルコード
そこで，ループの一部を分割し，最内ループのループ長のうち最大になるものを検索する．次に，最内ループの終端を検索した最大ループ長に変更
し，ループ間の依存関係を解消する．最大ループ長に満たないループはマスクにより処理を行わないようにしている．依存関係が解消したため，ループ交換を行い，外側のループを最内側ループに移動することでループ長を長く確保し，ベクトル化率を高めた(図 5)．さらに，ベクトルプロセッサに搭載されている ADB を有効に活用できるように最内ループのブロック化を行う． SX-ACE は１

図 6 サブルーチン C の高速化後のコード命令で 256 要素を同時に演算できるため，ブロッ
ク長は 256 とした(図 6）．

高速化後の性能を表 4 に示す．高速化によりベクトル化率が 12.23%から 98.91%に改善し，ベクトル長が 12.3 から 151.3 に伸びることで，サブルーチン単体で約 6 倍の性能改善を得ることができた．

表4

サブルーチン C の高速化の効果

実行時間 ベクトル ベクトル 加速率
化率
長
[sec]
before 1164.311
13.23
12.3
6.11
after
190.607
98.91
151.3

以上の高速化を適用後のコード全体の性能を表5 に示す．シングルコア実行で高速化前に比較し約 17 倍の性能改善を得ることができた．

表5before
after

高速化後の性能

実行時間[sec]
13397.172
790.884

加速率16.94

図 7 並列化性能
4.4 並列化
本コードはベクトル化を行うことでシングルコ5 まとめアでの性能を高めてきた．さらなる性能の向上の

本稿では反応・相変化を伴う多分散系混相流シ

ためには，ネットワークに接続される複数のノーミュレーションコードに，本センターが今までにドを利用し，多数のコアで計算をする並列化を検取り組んできた高速化支援で得られた高速化技法討する必要がある．並列化にはノード内のコア毎を適用し高速化を行った事例について紹介した．に プ ロ セ ス を 処 理 す る MPI(Message Passing

本センターの提供しているスーパーコンピュータInterface)並列を行う Flat MPI を用いた．

SX-ACE の性能を最大限に引き出すためには，ベク

本コードの計算領域は 1,901×161 の 2 次元空間トル化率を高く保つ必要がある．また，SX-ACE はで，この空間中に粒子が存在する．この 2 次元空従前のベクトル型スーパーコンピュータより，さ間中の熱や圧力等の計算後，粒子の相互作用の計らにノードを多数接続することで性能の向上を図算を行う．並列化をするにあたり，計算領域を分る構成となっている．従って，今後の高速化支援

割し分割後の各領域に MPI のプロセスを割り当て，においては，シングルコアでの高速化に加え，複各領域の熱や圧力等と存在する粒子に対して相互数のノードのコアを利用する並列化を行うことで作用の計算を行うように並列化を行った．

性能を引き出していくことが重要となる．

図 7 に並列化後の性能を示す．2 次元分割を行

本コードは，高速化により数時間程度を要する

うと計算の初期で粒子が存在しない領域が多数存

シミュレーションを数分程度で計算することが可在し，プロセス間のインバランスが大きくなるた

能となった．今後，さらなる大規模なモデルでシめ，計算領域を x 方向(1,901)での 1 次元分割とすミュレーションを行うことで安全防護技術の発展

る．高速化後のシングルコア実行に比較し，64 プに貢献していきたい．

ロセス実行で約 9 倍の性能改善を得るにとどまった．実行後，時間が経過すると，各領域に存在する粒子数に偏りが発生し，プロセス間のインバランスにより，性能向上が抑制されていると考えらえる．今後，さらなる高速化のため，多数の粒子が存在する領域をさらに細かく分割することや，ノード内の並列化に OpenMP を用いた Hybrid MPIなどの検討が必要である．

参考文献[1] 東北大学サイバーサイエンスセンター， 高速化推進研究活動報告，第 6 号，p8-12，(2015)．