OpenFOAM の性能強化と 1 千億格子規模データのポスト処理の試み 清水建設株式会社 PHAM VAN PHUC 内山学
京 における OpenFOAM に関する取組み 第 1 回 OpenFOAM ワークショップ (2013) 京 への移植 10 億格子計算の壁 解決策 (Int32, プリ ポスト ) 第 2 回 OpenFOAM ワークショップ (2014) 1 万並列計算の壁 解決策 (MPI プラットフォーム ) 第 3 回 OpenFOAM ワークショップ (2015) 超並列 超大規模解析の性能評価 (10 5 並列 1 千億格子 ) 第 4 回 OpenFOAM ワークショップ ( 本日 ) OpenFOAM 性能強化, 超大規模ポスト処理 (1 千億格子 ) 2
内容 OpenFOAM の性能強化 ( 内山 ) 1 千億格子規模データのポスト処理の試み ( フック ) 3
OpenFOAM の性能強化ー thread 並列化 ( 非構造格子対応 ) - 清水建設株式会社 内山学 4
1. 目的と背景 目的 : 流体解析やFEM 解析をHybrid 並列化する ( 本報告ではthread 並列化まで ) 背景 Flat MPIによる超並列計算では通信の負担が大 流体コードで, 直交格子については文献 [3] で検討したが 実用面では非構造格子に適用できる方法が必要 格子単位のマルチカラー法では PCG 法の収斂性が悪い 色数を増やせば収斂性は改善されるが, 並列性は悪化 係数行列の行単位の並列化は性能が悪い( 演算量小 ) 使用コード :OpenFOAM-v1606+(pisoFOAM) [3] 内山, 他, OpenFOAM による流体コードの Hybrid 並列化の評価, 第 151 回ハイパフォーマンスコンピューティング報告発表会,2015. 5
2. 格子のオーダリングと行列格納方法 格子のグループ分け Figure 1 32 32 格子のグループ分け 6
2. 格子のオーダリングと行列格納方法 モデル :motorbike グループ数 :(32, 4, 2) Figure 2 係数行列の非零項の分布イメージ 格子の 99% 以上が level=0 と level=1 に含まれる 7
2. 格子のオーダリングと行列格納方法 対角ブロックのオーダリング U[*][4] の形に格納 (m1=4 の場合 ) Figure 3 対角ブロック内のオーダリング 8
3. 計算方法 対角ブロックの A 部分 (m1=3): 前進代入 A 部分の行数 対角ブロックの A 部分 (m1=3): 後退代入と行列ベクトル積 9
4. 数値実験 : 高層ビル 風速 10m/s 1 Δt=0.0001 s 最初の 2 step の圧力方程式を解く ICCG 法を計測 計算機 :Xeon E5-2697 v2 (2.70GHz, 12 cores, 30 MB Cache) 2 10
4. 数値実験 : 高層ビル 2 (6 色 ) CM: Cuthill-Mckey ordering S: profile reduction by Sloan EGMC: multi-color b128_s0: level0 が 128 ブロックブロック内が Sloan b128_s2: b128_s0 に対してブロック内を部分的に U[*][m1] と保持 number of iterations 4000 3500 3000 2500 2000 1500 1000 CM0 EGMC b128_cm0/1 b128_s0/2 b256_s2 500 0 11
4. 数値実験 : 高層ビル 収斂計算回数の増加以上に遅い 計算の制御が増加並列時に逆転 16 threads での高速化率が良くない 対角ブロック内 Sloan で悪化 対角ブロック内提案方法で高速化 Figure 6 Elapsed Time of ICCG Figure 7 高速化率 12
4. 数値実験 :motorbike 風速 10m/s 3 Δt=0.0001 s 最初の 10 step の圧力方程式を解く ICCG 法を計測 13
4. 数値実験 :motorbike 4 CM: Cuthill-Mckey ordering S: profile reduction by Sloan b128_s0: level0 が 128 ブロックブロック内が Sloan b128_s2: b128_s0 に対してブロック内を部分的に U[*][m1] と保持 16 threads での高速化率が良くない Figure 9 Elapsed Time of ICCG Figure 10 高速化率 14
4. 数値実験 EGMC ICCG 法の収斂性は著しく悪化し, 高速化率も余り良くない b128_s2 CM0と比べてICCG 法の収斂性を殆ど悪化させない Thread 並列化しない状態で計算時間を約 5% 短縮 8 threads 時で高速化率が5.3であり 良い結果と言える 16 threadsでは高速化率の伸びが良くない 15
他の部分の高速化の例 surfacescalarfield phihbya ("phihbya", (fvc::interpolate(hbya) & mesh.sf()) + fvc::ddtphicorr(rau, U, phi)); コードを展開 : 多次元配列化 作業配列削除 #pragma omp parallel if(nths>1) { #pragma omp for for (register idata i=0; i<ncells; i++) { irauu[i][0] = irau[i]*(c0*iu0[i][0]-c00*iu00[i][0]); irauu[i][1] = irau[i]*(c0*iu0[i][1]-c00*iu00[i][1]); irauu[i][2] = irau[i]*(c0*iu0[i][2]-c00*iu00[i][2]); } #pragma omp for nowait for (register idata fi=0; fi<nfaces; fi++) { idata i=p[fi], j=n[fi]; rdata s1 = isf[fi][0]*( iweights[fi]*(iu0[i][0]-iu0[j][0]) +iu0[j][0] ) +isf[fi][1]*( iweights[fi]*(iu0[i][1]-iu0[j][1]) +iu0[j][1] ) +isf[fi][2]*( iweights[fi]*(iu0[i][2]-iu0[j][2]) +iu0[j][2] ); rdata ddtcouplingcoeff = 1.0 -min( fabs(iphi0[fi]-s1)/(fabs(iphi0[fi])+vsmall), 1.0 ); 当該箇所の計算時間 : 高層ビル (9,973,802 格子 ) 2.21 sec 0.66 sec ( コードを展開 ) 0.35 sec (8 threads) 未展開部分があり十分に並列化されていない こういう部分はたくさんある C++ の良さはなくなるが速くなる } 以下略. rdata gamma_ra = iweights[fi]*(irau[i] -irau[j]) + irau[j]; rdata s2 = isf[fi][0]*( iweights[fi]*(irauu[i][0]-irauu[j][0]) +irauu[j][0] ) +isf[fi][1]*( iweights[fi]*(irauu[i][1]-irauu[j][1]) +irauu[j][1] ) +isf[fi][2]*( iweights[fi]*(irauu[i][2]-irauu[j][2]) +irauu[j][2] ); rdata ddtphicorr = rdt *ddtcouplingcoeff *( gamma_ra*(c0*iphi0[fi]-c00*iphi00[fi]) -s2 ); rdata s3 = isf[fi][0]*( iweights[fi]*(ihbya[i][0]-ihbya[j][0]) +ihbya[j][0] ) +isf[fi][1]*( iweights[fi]*(ihbya[i][1]-ihbya[j][1]) +ihbya[j][1] ) +isf[fi][2]*( iweights[fi]*(ihbya[i][2]-ihbya[j][2]) +ihbya[j][2] ); iphihbya[fi] = s3 + ddtphicorr; 計算機 :Xeon E5-2637 v3 (3.50GHz, 4 cores) 2 16
1 千億格子規模データのポスト処理の試み 清水建設株式会社 PHAM VAN PHUC
超並列 超大規模解析の性能 第 3 回 OpenFOAM ワークショップ (2015) Weak Scaling ( 単相流解析ソルバ pisofoam) PHAM VAN PHUC 内山学 ( 清水建設 ), 井上義昭 浅見暁 (RIST), 千葉修一 ( 富士通 ) : 京 コンピュータでの C++ 型流体コードにおける MPI の評価, 第 151 回ハイパフォーマンスコンピューティング研究発表会論文集 (2015.9)
超並列 超大規模解析の性能 第 3 回 OpenFOAM ワークショップ (2015) Weak Scaling ( 単相流解析ソルバ pisofoam) 1,000 億格子規模の計算は可能である
最新 OpenFOAM の並列性能 Weak Scaling ( 二相流解析ソルバ interfoam) 清水建設 理研 富士通との共同研究 計算コア部分の大きな技術課題はほぼクリア プリ ポストの課題 32 億格子 256 億格子 PHAM VAN PHUC( 清水建設 ): Large Scale Transient CFD Simulations for Buildings using OpenFOAM on a World s Top-class Supercomputer, The 4th Annual OpenFOAM User Conference 2016/ 2016 North American OpenFOAM User Conference, Keynote 講演, 2016
プリ ポスト処理の課題 シリアル処理の限界 ( データハンドリング メモリ不足 ) データの分割 結合が非常にかかる 10 億格子の壁 10 億格子データ 1TB メモリの利用 京のプリ ポスト処理 PC: 1TB モデル作成 ( プリ処理 ) 初期化領域分割 シミュレーション データ結合ポスト処理 出力 可視化 Computational domain in a single processor Computational domain, results in a single processor
プリ ポスト処理プロセスの提案 第 1 回 OpenFOAM ワークショップ (2013) 分散処理による解決 OpenFOAMの既存ツールの改良モデルの作成 ( プリ ) 初期化ロードバランス シミュレーション ポスト処理 データ出力 可視化 画像処理 データ処理等
京 における可視化プロセスの実装 風洞実験の再現 (64 億格子 ) 広域市街地の再現 (20 億格子 ) 1,000 億格子規模データ可視化のトライ
大規模 超高速の可視化システムの開発取組み TSUBAME2.0 スパコンでの実施検討 (2012) PHAM VAN PHUC と渡辺宏一 ( 清水建設 ): 建築環境の超高速計算可視化システム, 招待講演, 2012 http://www.cybernet.co.jp/avs/seminar_event/conf/18/program.html#03 24
シミュレーション技術のコンセプト 精緻なデータの利用 ( プリ ) 3 次元データの利用 高解像度計算モデルの作成 プロ 素人とは関係なく 自動作成 高速 大規模計算の実施 ( 解析 ) 最新ハード 分散処理 超高速 超並列計算の実施 高精度予測 リアルタイム性処理 どこでも活用できる ( ポスト ) リモート クラウド技術の活用 客先と対話可能
シミュレーション技術のフレームワーク Database (3D データ ) どこでも活用できる ユーザー Platform 画像レンダリング 流体解析 構造解析 IT Services 京 : 解析の高精度化 CPUs GPUs CPU GPU クラウド ユーザー
数値流体解析データ容量の次元 2K 4K 8K 1KB~100MB 1GB~1TB 1TB~100TB 画像レンダリング データ合体 Screen pixels Image pixels 必要な可視化データ 計算領域物理情報 可視化エンジン Paraview Ensight FieldView
戦略的な可視化方法の概要 2K 4K 8K 1KB~100MB 1GB~1TB 1TB~100TB 画像レンダリング 低解像度化データ伝送 Screen pixels 画像伝送 Image pixels 必要な可視化データ データ合体 画像処理 計算領域物理情報 Multi-Screen pixels 多並列の画像レンダリング LOCAL (PC) REMOTE ( 京 コンピュータ等 )
京 での可視化エンジンの移植 対応ソフト ParaView-3.12.0 ParaView-3.98.1 ParaView-4.1.0 コードの一部修正 移植 画像処理ソフトは独自コード
対象とする解析モデル (1,000 千億 ) 京 利用課題の 竜巻シミュレーション 2m 竜巻のコア 建物 16m processors: 256 24 16=98,304 総格子数 : 25,856 2,424 1,616 1 千億 3m
データ容量 ハンドリング 計算格子 : 3.19TB *) 物理量 (U, p) : 0.56TB/ ステップ *) データファイルのハンドリング Processor 数 : 98,304 ステージングイン時間 : 約 10 分 (proc=98,304) *) tgz 圧縮フォーマット
画像の概要 Resolution 解像度 X Y Pixels per Frame ( 画素数 ) Full HD(2K) 1,920 1,080 2,073,600 4K 3,840 2,160 8,294,400 8K 7,680 4,320 33,177,600 16K 15,360 8,640 132,710,400(1 億 ) 32K 30,720 17,280 530,841,600(5 億 ) 48K 46,080 25,920 1,194,393,600(11 億 ) 64K 61,400 34,560 2,123,366,400(21 億 ) 世界最高画素数カメラ : 5,060 万画素 ( キャノン EOS 5Ds) 超高精細ディスプレイ : ~8K 世界最大画像サイズ : 3,650 億画素 ( 画像処理 )
可視化の画素数 ( 投影面 ) 1pixel/ 格子 4pixels/ 格子 9pixels/ 格子 十分な画素数を用いる必要がある
全メッシュ図の可視化 9pixels/ 格子 フォト表示サイズ :722.5KB Pixel サイズ :91,648 17,184 ~ 16 億画素数 ~ 64K 34
建物周辺のメッシュ図 9pixels/格子 サイズ:2.21MB Pixelサイズ 8,000 8,000 ~ 6400万画素数 8K 35
風速コンター図の可視化 9pixels/ 格子 フォト表示サイズ :16.3MB Pixel サイズ :91,648 17,184 ~ 16 億画素数 ~ 64K 1,000 億のデータを技術的に可視化できることを確認! 36
計算格子の低解像度化による変化 1,000 億格子 データマッピング 2 億格子 画像はやや劣化しているが 現象の観察について十分
まとめ 1,000 億格子規模データのハンドリングや可視化は分散処理により行うことが可能である 現象の把握を目的とした場合 部分可視化やデータの低解像度化を行うことは有効な方法である 可視化レベルによってそれを限らない