RIST ニュース No.64(2018) GPU の気象気候海洋モデルへの適用とその最適化を目指した取り組み ActivitiestowardsanapplicationofGPU tonumerical weather/climate/oceanmodelanditsoptimization

Size: px

Start display at page:

Download "RIST ニュース No.64(2018) GPU の気象気候海洋モデルへの適用とその最適化を目指した取り組み ActivitiestowardsanapplicationofGPU tonumerical weather/climate/oceanmodelanditsoptimization"

はすなみおか
5 years ago
Views:

1 GPU の気象気候海洋モデルへの適用とその最適化を目指した取り組み ActivitiestowardsanapplicationofGPU tonumerical weather/climate/oceanmodelanditsoptimization 高度情報科学技術研究機構山岸孝輝次世代のエクサスケールスーパーコンピュータでは性能及び電力効率の高さが求められるその両方を備えたシステムとして演算性能とメモリ転送性能に特化したアクセラレータを用いたヘテロジニアスなシステムが広がりを見せている演算性能電力効率の世界ランキングの両方においてアクセラレータを用いたシステムが多数上位に入っておりその性能に対する期待は大きいその一方アクセラレータは数値計算に特化したプロセッサであることに加え CPU との複合システムであることからユーザの開発及び最適化の負担が大きい事が問題となっている本稿ではこれらのアクセラレータの内最も普及しているNVIDIA 社のGPU を取り上げ現在及び今後のスーパーコンピューティングにおけるGPU 適用の背景基礎並びに今後の課題について概説し併せて我々が行っている気象気候海洋モデルへのGPU の適用とその最適化について紹介するはじめに地球温暖化を初めとする気候変動の予測極端現象異常気象による災害の予防と即時対応の手段としてまたそれらの基礎となる物理プロセスを理解する手段として気象気候及び海洋の数値シミュレーションの重要度は高まる一方であるそれ故気象気候海洋モデルの高度化の需要は大きい本稿では気象気候海洋モデルの有効な高度化手段の一つとしてアクセラレータ (NVIDIA GPU) への対応についてその背景基礎並びに今後の課題について概説し併せて我々の研究成果を紹介する 2. 気象気候海洋モデルの特徴と開発の背景気象気候海洋モデルでは複雑かつ多くの種類の現象をまとめて扱う事が特徴である例えば気象の数値モデルで扱う大気現象をスケール ( 代表的な大きさ ) の違いで考えると空間スケールは地球 1 周 4 万 km 大気拡散の運動は数 mmである季節の移り変わりは数ヶ月気候変動となると数十年風の流れは数時間から数分であるつまり時間空間スケールの幅は非常に大きい関連する物理現象も大気放射エアロゾル ( 大気中の塵 ) 大規模凝結積雲対流降水に地表面過程など多岐にわたりこれらが相互作用しながらいわゆる気象システムを形成している気象の数値モデルはこれらの現象の相互作用を時間空間的に広いスケールで扱う高度な数値モデルといえるより良い精度の予報並びに予測とそのプロセスの理解にはモデルをより高解像度で実行することが求められる加えて防災という観点からは予報の即時性が求められるため予報が得られるまでの経過時間の減少も重要な要素であるまた統計的手法により予報の精度を上げるためには十分なアンサンブル数の確保が必要であり経過時間の減少

2 -30- は重要である数値モデルはその名の通り実際の自然現象を近似して数値化しているのであるがより詳しく現実を再現理解するためにはなるべく近似を行なわずより多くの力学物理過程を直接表現することつまりモデルの精微化が必要となる以上に挙げた高度な気象気候海洋モデルへの要求事項 ( 高解像度化高速化精微化 ) に対する制約の一つは計算機の能力と実効性能である例えばメッシュサイズに依って数値モデルに含む要素やモデル化の手法も決まりその結果研究の対象にも制約が生じる計算機の能力実効性能と研究対象の関係は線形では無くある閾値を境に質的な変化をも含むこれは同じ計算機の能力でも高度化の内容次第で研究の質ががらりと変わりうるということである気象気候海洋モデルは社会的な重要度が大きくかつ影響を与える分野が広いため利用者の数も使い方も自ずと増える直接開発に携わらない人や専門的な領域のみ詳しい人が占める割合が増える事に加え数値モデルを実行する環境も多様化するそれ故なるべく汎用性が高くて実行環境を選ばないこと簡単に開発高度化できることコードの理解に誰もが簡単に取り組めること等が重要となるこれまでに挙げた課題を計算科学の言葉で言い換えると気象気候海洋モデルの高度化は弱いスケーリング ( 高解像度化 ) と強いスケーリング ( 高速化 ) の両方が求められていてそのコードはさらなる複雑化を求められ生産性 ( 容易な開発 ) やポータビリティ ( 多様な実行環境に対応 ) も求められるという実に挑戦的な取り組みといえるこれらの課題に対してどのように対応してきたかについてこれまでの日本での気象気候海洋モデル開発の歴史を振り返ってみるとベクトル計算機の並列クラスタである地球シミュレータから超並列なスカラ計算機への移行を経て京コンピュータでの大規模実行を可能とする開発と実行など革新的かつ多大な努力が多方面から払われてきたこの開発の流れに今後さらに影響を与えるものに電力効率改善の要求があるその要求に対応するため次世代のエクサスケールスーパーコンピュータではこれまでの超並列化の流れに加えてアーキテクチャやシステム構成の変更をも伴う大きな転換を求められている電力効率を上げる方法としてプロセッサの動作周波数を下げると共にコア数を増やし ( 半導体プロセスの集積度の増加分をコア数増加に転換させて ) 電力消費量を抑えたままでチップあたりの演算性能は向上させる方法と目的に特化した設計として無駄を省く方法があるベクトル計算機からマルチコアを搭載した超並列なスカラ計算機への移行においては前者の方法が主に取られてきた今後の流れではさらに後者の方法も重要とされている超並列化及び構造の特化で無駄を省く方法の両方を採用するプロセッサとして近年広がりを見せているのは従来のCPU と組み合わせて用いるアクセラレータである CPU とアクセラレータの組み合わせを採用したヘテロジニアスな構成となるスーパコンピュータ ( 以下スパコン ) は演算性能 (TOP500) 電力効率(Green500) の世界ランキングの両方において多数上位に入っており TOP500 での最新の1 位はGPU を搭載した米国のSummit である日本では産業技術総合研究所のABCI(5 位 ) 東京工業大学の TSUBAME3.0(19 位 ) などがランク入りしているアクセラレータの広まりと同時に課題として上がってきたのはプロセッサの超並列化とスパコンシステムのヘテロ構成化 (CPU+ アクセラレータの組み合わせ ) がアプリ制作側の負担となることであるこれは気象気候海洋モデルのみが抱える問題では無い

3 本稿では特に気象気候海洋モデルがこれから迎える超並列化の時代に抱えるであろう課題とその対応について関連する背景や基礎を紹介しつつ纏めていく後述するが気象モデルで扱うカーネルは構造格子でメモリアクセスも系統的である故に他のアプリにも適用可能な高度化のノウハウを含むものであり他の事例にも適用可能と考えているまた近年最も広まっているアクセラレータはNVIDIA 社のHPC 向けグラフィックカード (GPU) を採用したものであり本稿では全てNVIDIA 社のGPU を前提とする 3.GPU とは本章では GPU の基礎について気象気候海洋モデルに関係が深い部分を中心にまとめるさらなる解説等はNVIDIA の公開文書 [1,2] を参照されたい元来はグラフィック用のプロセッサを HPC 向けに利用したもので HPC で主に使われる上位機種ではグラフィック処理ではさほど必要としない倍精度演算が強化されており浮動小数点の表現ではIEEE に準拠しデバイスメモリにはECC が搭載されているなど科学技術計算において求められる機能を備えているコンパイラやライブラリプロファイラ等のソフトウェアも提供されハードソフトの両面で十分な機能を備えているこれまでの一般的なCPU に比べて数値計算に特化したことで性能に対する電力の効率が高いことが特徴であるしかしながら特化したが故にCPU とは異なる並列化手法が求められプログラミングにおいていくつかの制約が発生する数値計算に特化して構造を単純化させたことは使いにくさと使いやすさの両方につながる GPU は基本的なやり方を守って実装しないと十分な性能を出すことができず理論演算性能で劣るCPU よりも実効性能が出ない事が殆どであるこれは多くのユーザが抱える問題ではあるが筆者は単純であるが故の利点例えば使いやすさもあると考える構造が単純である場合ハードの機能に依る余計なレイテンシが発生する機会が少なくなり ( 例えばGPU のL1 キャッシュはコヒーレンシを取らない ) 加えてモデルの振る舞いを把握しやすく書き手が操作しやすくなるモデルの最適化では性能の阻害要因を特定しその原因を調査した上でコードを修正していく複雑な構造の場合原因を特定するためのコストが大きくなり加えて書き手の狙い通りに修正することが難しいこの欠点と利点のどちらが大きいかはアプリの特性や書き手のGPU に対する習熟度によるが筆者らは気象気候海洋モデルについては先行研究や筆者らの事例では良い結果が出ていることから利点が上回ると考える立場であるここで GPU の実行モデルについて要点を纏めるデバイスであるGPU はホストであるCPU から独立しており GPU メモリの確保と転送はCPU から行う GPU 上のメモリをCPU から確保した上で CPU メモリからGPU メモリにデータを転送する ( ホストデバイス転送 ) GPU カーネルは転送されてきたデータを更新した後 CPU メモリに転送してGPU での処理を終える GPU カーネルでの処理は CPU ではデータが連続する方向にループを回転させることで多数のデータを連続して処理していくが GPU は1つのデータを1つのスレッドに割り当て CPU でのループの回転をスレッドの割り当てに置き換えるこのように多数のスレッドをGPU が持つ数千の演算コアにデータ並列性を活かすように割り当てていくことで多数のデータの処理を行う ( 図 1) -31-

図 1 GPU の実行モデル科学技術計算にてオーバヘッドまたはメモリの読み込みなど待ちが発生することは大きな問題であるが CPU はこの問題をキャッシュとその他の機能 ( 例 :SIMD 処理プリフェッチ分岐予測 ) で対応する一方 GPU は多数のスレッドを立て待ちの発生時にはスレッドを切り替えることで対応する GPU

4 図 1 GPU の実行モデル科学技術計算にてオーバヘッドまたはメモリの読み込みなど待ちが発生することは大きな問題であるが CPU はこの問題をキャッシュとその他の機能 ( 例 :SIMD 処理プリフェッチ分岐予測 ) で対応する一方 GPU は多数のスレッドを立て待ちの発生時にはスレッドを切り替えることで対応する GPU は処理の粒度を極力細かくした大量のデータを扱うことを前提しておりあるスレッドが待ちの時は実行の準備が完了した他のスレッドに切り替えて先に実行させてしまうその実行の裏で待ちの処理をさせて次のスレッドの切り替えに備えさせる裏を返すと GPU は高いデータ並列性が無いと待ちのレイテンシを隠蔽することができず無駄なコストの増加につながるこのように GPU とCPU は実行モデルが全く違うものであるシリアル処理が多い問題に対しては CPU が得意とするそれに対しデータ並列性が高い処理を含む問題に対しては GPU が得意とする GPU 実装で問題となるのは CPU とGPU の複合システムであるが故に物理メモリが 2つ存在しその間で転送が必要なことである転送を記述することでプログラミングが煩雑になりレイテンシを含む転送コストはパフォーマンスに直接影響する現状 GPU に対応する気象気候海洋モデルは初期処理をCPU で行いメインループ内の演算処理は全てGPU で実施しているまた MPI 通信やファイル入出力処理は本質的に CPU のみでしかできないためホストデバイ -32- ス転送は避けられないこれらに対しては GPUDirectRDMA 通信ストリームを用いた非同期処理などの手法でコストの増加を抑える事が可能である ( 手法の詳細は [3] を参照 ) 現状の実装例の殆どでは気象気候海洋モデルのメインループが回転している間は MPI 通信やファイル入出力処理以外ではCPU はアイドル状態にあり見方を変えるとCPU 資源の無駄が生じているといえるアプリケーションをGPU 対応させる方針として 3つ挙げられる GPU 対応ライブラリへの置き換え指示行挿入またはCUDA による直接的な記述であるライブラリの利用では通常のCPU コードで既存のライブラリを利用する場合とほぼ同様に利用できるこの場合ライブラリのコール時にホストデバイス間の転送を伴うためその転送コストとGPU による高速化のどちらの効果が大きいのか正しく評価する必要がある指示行挿入によるものでは仕様の更新改良の頻度が高く最も広く使われているものにOpenACC がある [4] OpenACC は OpenMPと同じく並列化のオープンな規格であり複数のメーカによってそれを実装したコンパイラがリリースされている機能豊富でライブラリ利用よりも実装での自由度が高い利用方法はOpenMPと同様並列化するループの前後に指示行を挿入することでそのループをGPU での並列化に対応させる OpenACC 対応コンパイラのオプションにて指示行の有効無効を指定できるこれにより CPU 用 GPU 用の実行形式を同一のソースコードから作成することができるためコードの可搬性が高い実装方法といえる後述するCUDA に比べて容易に実装可能であり開発メンテナンスのコストを削減できるしかしながらオープンな規格故にNVIDIA GPU に特化した機能を含んでおらずハー

5 ドが持つ全ての機能を活用できない例えばシェアードメモリの活用など後述する CUDA にてその有効性が示されている手段を取れないという問題点がある最も自由度が高い実装方法に CUDA を用いた記述がある CUDA はNVIDIAGPU 向けの言語拡張で現在主なものではC 言語向けとFortran 向けがリリースされている GPU のスレッドを陽に扱う事ができるため柔軟なコーディングが可能となりより深い最適化が可能になるしかしながら CPU 版のコードからポーティングする場合コードの大半を書き換えることが必要となるため CPU GPU 間でのコードの互換性を保つことは難しいこのように実装の容易さと性能の出しやすさは相反するものとなっているこれら3 つの実装方法はそれぞれ異なる特徴を持つが実装において排他的では無く併用した実装が可能であるアプリに応じて選択し適切に当てはめていくことが求められる 4. 気象気候海洋モデルの基礎 4.1 モデルの構造とモデル化気象気候海洋モデルの基本的な構成モデル化の方法数値解法などの基礎のうち GPU 向け実装に関係が深い部分を中心にまとめる気象気候海洋モデルは大きく分けて力学過程と物理過程の2つから成る気象モデルの場合力学過程は流体の運動方程式を基本として連続の式熱力学の式水蒸気の式などで構成される物理過程は力学以外の要素の総称で気象モデルの場合大規模凝結積雲対流大気放射計算エアロゾルや地表面過程などが含まれる以上から気象気候海洋モデルは非常に複雑な非線形数値モデルといえる力学過程は水平方向の差分化でスペクトル法と格子点法に分類されるが隣接ノード間での通信が主となる格子点法を採用する気象モデルが増加しており本稿で紹介する気象モデルは全て格子点法を用いている格子点法を用いて力学過程を構成する方程式を離散化するにあたり高次の離散スキームを取ることは少なく隣接する格子点との間で離散化されるスキームとなることが殆どである構造格子かつ系統的なメモリアクセスではあるものの隣接格子への3 次元方向での参照が多く加えて物理過程に比べて複雑な演算を含まないことからメモリバンド幅ネックになりやすい物理過程は各格子で計算が閉じるか鉛直方向のみで依存するものが大半である力学過程に比べて処理が複雑となり加えて特殊関数を多く含むなど処理の中で演算が占める割合が多く演算ネックになりやすい以上から気象気候海洋モデルは計算特性が異なるものを含む計算科学の視点からも複雑な数値モデルといえる自然現象全てを微細なプロセスまで全て計算しているわけでは無くより大きなスケールの変数から経験的に微細なプロセスの役割を見積もるいわゆるパラメタリゼーション化がなされている例えば積雲が発生する水平の空間スケールは数百メートルであるため数十キロの格子間隔の気象モデルでは表現することができないその場合その数十キロの格子を代表する物理量のうち積雲の発生と関係が深いものから物理的経験的な手法によってその格子にて発生する積雲の平均的な値を見積もることで小さいスケール ( 積雲 ) が大きいスケール ( 数十キロの格子の状態 ) に与える影響を評価しているこの手法は積雲パラメタリゼーション化と呼ばれているパラメタリゼーション化の存在は予報にある程度の不確定性を与え得るため解決策として計算の空間解像度を上げて直接計算を行いパラメタリゼーションを外してしまうことが挙げられる水平解像度を細かくして積雲の発生を直接計算する全球非静力学 -33-

6 -34- 気象モデルの開発も行われてきた [5] 直接計算するか否かでモデルの結果には質的な違いが生じそれまで表現できなかった物理を捉えることも多く計算能力の向上効率化がもたらす気象気候海洋モデルの精微化に対する需要は大きい気象気候海洋モデルは以上に挙げた全ての要素が相互作用しながら時間発展していく非線形問題である初期データ ( 予報変数の初期値境界条件 ) を読み込んだ後は全ての要素を含んだ時間ループを何度も繰り返していくしかしながら全ての要素に適切な予報変数を設定しそれらの連立方程式を解くのは計算量が膨大となってしまう実際のモデル実装ではいくつかの仮定の下で近似が行われ時間スケールや処理の性質に応じて分割された上で計算される例えば力学過程 ( 風速地表面気圧雲水など ) を計算する際物理過程 ( 大気放射エアロゾル積雲対流降水や地表面過程など ) は更新せずに一定としているそして更新した力学の場を元に物理過程のそれぞれの要素を逐次に更新する各要素の更新の順番組み合わせ更新の頻度をどのように設定するかは気象気候海洋モデルにとって不確定要素の1つといえるそれは予報の精度並びに現象の解釈に影響を与えるのみならずモデルの実行速度にも関係しているそれ故気象学計算科学の両方の視点から気象気候海洋モデルの構成を検討する必要がある 4.2 気象気候海洋モデルへのGPU 適用事例求められる即時性という社会的な重要性と構造格子でメモリアクセスが系統的という数値モデルとしての性質の良さから GPU に対応させる取り組みは早くから行われてきたしかし初期の事例では CUDA による部分的なポーティングに留まるものであったメインの時間ループに含まれる処理全てをGPU 上で初めて実行したのは CUDAC によるGPU 実装を行った気象モデルASUCA のTSUBAME2.0 での大規模実行の事例である [6] GPU 単体での性能改善やマルチGPU での通信の隠蔽等により高い実効効率の元での大規模実行を達成した CUDAC は最適化の自由度が高く高い性能を引き出すことができたが CPU 版コードとの互換性は保たれておらず開発メンテナンスにおける効率が良いものでは無かった文献 [6] で示した事例を含む初期の事例では CUDA を用いたGPU 実装で高い性能を示したものの研究室または研究室間での試みにとどまり気象予報等での実運用を前提とした開発ではなかった近年の開発ではコードのポータビリティ開発メンテナンスの工数削減の必要性から CUDA での開発から指示行挿入によるOpenACC を中心とした開発が主流となってきたまた大学研究機関企業での共同体制が多くのプロジェクトで確立されているこれは気象モデルのGPU 対応に十分な意義があることの現れといえる OpenACC を用いて大規模システムで多数ノードまで性能をスケールさせた例としては TSUBAME2.5 での NICAM 力学コアの OpenACC 実装 [7] の他スイス連邦工科大学スイス気象庁によるプロジェクト [8] にて気象モデルCOSMOを OpenACC での実装に加えDSL やライブラリとの併用で力学過程物理過程全てをGPU に実装し数千 GPU まで性能がスケールすることを示した [9] またスイス気象庁ではGPU スパコンであるPizKesch 上でGPU を用いた初の気象予報が行われており GPU が気象予報の実運用に耐えうることが初めて示されたこの他 NCAR 等の米国の研究機関が中心となって開発してきた領域気象モデルWRF や全球気象モデルMPAS では NVIDIA を含む

7 民間企業も参加した国際的なプロジェクトとして OpenACC を軸としたGPU への対応が進められている [10] 5. 計算科学としての気象気候海洋モデルの課題と今後ここまでは気象気候海洋モデルの GPU への適用と高度化についてその背景や基礎を概説した本章ではそれらを踏まえて課題と今後の対応についてまとめる気象気候海洋モデルでは特に力学過程においてメモリアクセスがカーネルの性能を律速する要因となることが多い一般に GPU ではメモリ転送性能の演算性能に対する割合がCPU に比べて低いため ( メモリ転送性能の絶対値ではCPU よりも優れているものの ) 演算の実行効率という視点からは問題となり得るまた CPU GPU 共に近年のスーパーコンピューティングではメモリ転送性能の演算性能に対する割合が低くなる傾向にあり気象モデルのこの性質は性能を向上させる上で大きな問題の一つとなっているまたカーネル内に一時変数を持つことも多くレジスタ消費量も多い特に物理過程に含まれる経験則に基づく式は複数のパラメータを含む事が多くレジスタ使用量が多くなりがちであるレジスタ消費量が大きいと同時に実行状態となるスレッドの数が減少しデータ並列性を確保できないまたレジスタスピルが発生してキャッシュ等への変数の一時退避が発生するなど性能を阻害する要因となるノード間通信は毎時間ステップにおいて複数回発生し通常のMPI で通信する場合はその都度 GPU 上のメモリからCPU 上のメモリに転送する必要があり CPU 単独での実行に比べて余計なコストとなり得るモデルの実行コスト分布に関してはコード全体で共通して使われるような処理が少ないため計算コストが集中するような箇所が存在しないことが特徴である全球スペクトル気象モデルの性能分析例では最もコストが集中したサブルーチンは全体の数 % 程度のコストでしか無く大半のサブルーチンは 1% 未満のコストしか占めなかったモデルの開発コストに関してはその行数が非常に多いこと複数の研究者による共同での開発であること並びにCPU での開発が基本であることが特徴である複数の要素から構成されるため必然的にコードの分量は大きくなるまた構成要素 ( 例えば積雲対流など ) それぞれを担当する研究者が独自の方針でCPU にて数値モデルを作成検証し持ち寄ったモデルを結合したものを1つの気象気候海洋モデルとして共有しているそれ故開発者全てがモデルの高度化に詳しいとは限らない気象気候海洋モデルのGPU での高度化は複数の課題を抱えておりその課題を解決するには計算科学の専門家のみでは不可能で地球科学を専門とする側との協力がいろいろなレベルで必要となるそのためには理解しやすく実装に時間がかからないプログラミング手法が必須といえるここまでコードの特性からのGPU 向け最適化の困難さ気象モデルコミュニティ故の開発に対する問題点の二つを述べたが後者についてはOpenACC が唯一の解である実際昨今はOpenACC での指示行挿入での開発がほぼ全ての気象気候海洋モデルで行われているしかしながら第 2 章で挙げた高解像度化高速化精微化全てに対応するには機能面の制約から限界があると言わざるを得ないまた OpenACC の元でGPU での実行を優先して最適化すると CPU 実行での性能が低下しかねない今後の取り組みとしてまず筆者が重要と考えることは OpenACC とCUDA を併用した実装である加えて将来的な案として先述した気象モデルの複合的な性質を元に気象 -35-

8 システムを構成する要素の並列化による高速化を紹介する OpenACC での高度化が難しい箇所は CUDA またはGPU 対応ライブラリを用いるその際 CPU 版コードと併存させることになりプリプロセッサなどの指定でコンパイル時に切り替えることが良いと考えるその部分については2 種類のコードが存在することになりメンテナンスのコストは2 倍となる気象気候海洋モデルは大量のコードを含むが長い間更新されないいわゆるレガシー化された部分もありそのような部分はCUDA 版とCPU 版のコードを併用して運用しやすい性能計測を行いコスト上位から OpenACC とCUDA を併用していくやり方もあるがコードの更新頻度も併用する部分を選択する指標の一つとなる先述の通り CUDA は最適化のスイートスポットが小さいが幸い気候気象海洋モデルは構造格子でかつ系統的なメモリアクセスであるため最適化の方針は立てやすい CUDA を用いた先行研究で高い実効効率を出したことがそれを示している OpenACC とCUDA を併用した実装の優先度が高いが将来的な案として気象システムを構成する要素の並列化による高速化の案を紹介する先述したとおり気象モデルは自然現象の内実際は同時に行われる要素 ( 例えば雲の凝結と大気の放射計算など ) をモデル化の段階で逐次処理に設定しているこの場合 GPU またはCPU で逐次に更新していくことになる複数の要素を並列に処理する様にモデル化することもその要素の性質によっては不可能では無いこれに基づき CPU 向けの処理とGPU 向けの処理を並行処理できるようなモデル化を検討する空いているリソース ( この場合はCPU) を活用することに加え CPU 向けの処理はCPU で高速に処理させることができればモデル全体の高速化につながるコアの単体性能が高くキャッシュの活用でレイテンシが大きい問題に対応可能な構造であるCPU と多数のコアによる多数の並列計算を可能にする GPU という異なる性格を併せ持つシステムは気象モデルのような複雑な要素を多く含むモデルには本質的に相性が良いといえる今後必要となるOpenACC とCUDA を併用した実装であるが OpenACC はCPU コード内のループをCUDA で記述されたGPU カーネルに置き換えるものであるため CUDA 実装の理解はCUDA での実装部分のみならず OpenACC での最適化にも役立つものであるそこで次章では非静力学海洋モデルを CUDA を用いてGPU 実装並びに高速化した筆者らの研究を紹介する 6.CUDAによる非静力学海洋モデルkinaco のGPU 実装と高速化非静力学海洋モデルkinaco のCUDA を用いたGPU 実装と高速化の事例を紹介する kianco の基本方程式等は [11] を参照されたい計算のメイン部分は前処理としてマルチグリッド法を用いたCG 法でPoisson 方程式を解く部分である高度化のさらなる詳細については参考文献 [12,13,14] を参照されたい CG 法前処理への混合精度の適用 [12] 気象モデルに含まれる浮動小数点演算は倍精度に設定されることが多い GPU は単精度演算のユニットと倍精度演算のユニットは別のハードで構成されていてチップ当たりの演算性能は単精度演算の方が優れているメモリ転送性能も単精度は倍精度に比べて実質 2 倍の性能に相当する高速に処理する上で倍精度から単精度への移行を検討することは十分な意味がある精度が落ちてしまうため計算結果への影響は評価する必要があり適用は慎重に行うべきであるここでは運動方程式のソルバである前処理付きCG 法の前処理のみに単精度での処理を適 -36-

用しソルバ本体であるCG 法は倍精度のままで実行したその結果収束回数は前処理含めて全て倍精度に設定時と同じ回数となり計算結果の精度も結果の解釈に影響を及ぼさなかったコアレスアクセスの促進 [12] Poisson 方程式に今回含まれる行列ベクトル積では各格子とその周囲 6

複数のスレッドが同時にアクセスするGPU ではストライド幅が 7つのアクセスとなってしまう ( 図 2) そこで配列要素の次元を並び替えコアレスアクセスとすることでメモリへのアクセス効率を向上させた図 2 コアレスアクセスの促進ードしレジスタに保存するシャッフル関数を用いて

それに対して連続にアクセスさせて別の物理量を計算していたこの計算方法は最初に計算した3 次元の物理量をキャッシュを介して次の連続アクセスで再利用することが目的であったが GPU の場合はレジスタを活用することでキャッシュよりもより高速にデータを再利用できる ( 図 4) 前述の2つシャッフル関数の活用

9 用しソルバ本体であるCG 法は倍精度のままで実行したその結果収束回数は前処理含めて全て倍精度に設定時と同じ回数となり計算結果の精度も結果の解釈に影響を及ぼさなかったコアレスアクセスの促進 [12] Poisson 方程式に今回含まれる行列ベクトル積では各格子とその周囲 6 点の参照を係数行列によって表現している合わせて7 点の係数は各格子固有の値をとるオリジナルのCPU コードでは各コア ( スレッド ) が逐次にメモリをアクセスしていくがその際に係数 7つがキャッシュライン上に連続して並ぶためキャッシュラインを有効に活用できるしかしながら複数のスレッドが同時にアクセスするGPU ではストライド幅が 7つのアクセスとなってしまう ( 図 2) そこで配列要素の次元を並び替えコアレスアクセスとすることでメモリへのアクセス効率を向上させた図 2 コアレスアクセスの促進ードしレジスタに保存するシャッフル関数を用いて他スレッドのレジスタに保存された値を自スレッドにロードして足し込み粗い格子を作成した図 3 粗い格子を作成するアルゴリズムループ融合とレジスタの活用によるデータ並列性の確保 [13] 物理過程に含まれるカーネルにおいて CPU コードではある3 次元の物理量を計算しそれに対して連続にアクセスさせて別の物理量を計算していたこの計算方法は最初に計算した3 次元の物理量をキャッシュを介して次の連続アクセスで再利用することが目的であったが GPU の場合はレジスタを活用することでキャッシュよりもより高速にデータを再利用できる ( 図 4) 前述の2つシャッフル関数の活用 [14] CG 法の前処理に用いたマルチグリッド法では格子を粗くする細かくする操作が複数回実施される GPU でのこの操作は各スレッドによる細かい格子からのロードが非コアレスアクセスとなってしまう ( 図 3) キャッシュの活用により 2 回目以降のアクセスでは最適化されるものの本件ではさらなる高速化を目指しウォープシャッフル関数を活用したまずコアレスアクセスでロ図 4 ループ融合とレジスタの活用によるデータ並列性の確保 -37-

の処理を融合しレジスタを介して計算を各スレッドで閉じる同時にループ間の暗黙の同期が削除されるため 2 次元スレッドが独立に動きデータ並列性を活用できるブロック形状変更によるキャッシュブロック化の促進 [14] kinaco の格子形状の仕様から Poisson 方程式の係数が鉛直方向にてほぼ一様 ( 海面付近海底付近を除く ) という特性がある CPU コードでは同一の値となる係数は

[32,1,16] はおよそ2 倍高速である (-1: n1+2,-1:n2+2,-1: n3+2) としてメモリを確保しているモデルの仕様上 n1 らは2の階乗となる場合が殆どでありその場合メモリへのアクセスでアライメントが適切なものにならずどの次元においてもトランザクションが最低で一つ無駄になる筆者らはオリジナルのCPU コードではMPI 通信時のオーバヘッド回避のために

10 の処理を融合しレジスタを介して計算を各スレッドで閉じる同時にループ間の暗黙の同期が削除されるため 2 次元スレッドが独立に動きデータ並列性を活用できるブロック形状変更によるキャッシュブロック化の促進 [14] kinaco の格子形状の仕様から Poisson 方程式の係数が鉛直方向にてほぼ一様 ( 海面付近海底付近を除く ) という特性がある CPU コードでは同一の値となる係数は同じメモリアドレスに格納してキャッシュの利用効率を上げる最適化がなされていた ( 図 5) この最適化はGPU のL1 キャッシュでも有効であるがキャッシュヒット率を上げるためにブロックの形状を鉛直方向の要素数が大きくなるように確保した同じブロックが処理を終了するまでストリーミングマルチプロセッサから解放されない特性を活用している評価例を挙げると [32,8,1] のブロック形状に対して [32,1,16] はおよそ2 倍高速である (-1: n1+2,-1:n2+2,-1: n3+2) としてメモリを確保しているモデルの仕様上 n1 らは2の階乗となる場合が殆どでありその場合メモリへのアクセスでアライメントが適切なものにならずどの次元においてもトランザクションが最低で一つ無駄になる筆者らはオリジナルのCPU コードではMPI 通信時のオーバヘッド回避のために転送データ ( 袖領域 ) を別途用意した1 次元配列に詰め替えた後に通信を行っていることに着目しこの1 次元配列をGPU カーネルから直接参照させて予報変数には袖領域を用意しないことにした ( 図 6) これによりメモリアクセスで無駄なトランザクション無くなることに加え受信した1 次元配列から予報変数への詰め替えのコストも削減出来た図 6 袖領域を持たない隣接通信図 5 ブロック形状変更によるキャッシュブロック化の促進メモリアライメント調整 ( ノード間通信用袖領域の削除 )[14] 各ノードでの予報変数は A(1:n1,1: n2,1:n3 ) としてn1 * n2 * n3 個の要素を有する場合袖領域は各次元で2 個つまりA コストが小さいカーネルの融合 [14] マルチグリッド法を用いた前処理では粗い格子 ( 例えば [8,8],[4,4] 程度 ) での演算が複数回実施されるその中に含まれる行列ベクトル積のカーネルは代表的なコストは経過時間にして数マイクロ秒程度であるがカーネルの起動に伴うオーバヘッドは 20 マイクロ秒程度とカーネル本体の計算に比べて無視できない大きさである反復法であるが故に呼び出し回数はモデル本体のステップ数の100 倍以上と非常に大きいそこでカーネルを融合することで起動に伴うオーバヘッドを削減した ( 図 7) オリジナルのCPU コードでは隣接間通信の前後に複数のカーネルが組まれていたが通信以外の -38-

カーネルを全て融合した余計なロードストアが減ったこととカーネル終了時にグリッド全体でとられていた同期がブロック内での同期となったことによるコスト減少の効果も含まれる図 7 コストが小さいカーネルの融合粒子追跡 [13] 算出した流速分布から海中の粒子 ( 栄養塩や放射性物質など ) 分布を決める粒子追跡計算にて CPU 版は粒子のデータ構造に連結リストを使用していた

11 カーネルを全て融合した余計なロードストアが減ったこととカーネル終了時にグリッド全体でとられていた同期がブロック内での同期となったことによるコスト減少の効果も含まれる図 7 コストが小さいカーネルの融合粒子追跡 [13] 算出した流速分布から海中の粒子 ( 栄養塩や放射性物質など ) 分布を決める粒子追跡計算にて CPU 版は粒子のデータ構造に連結リストを使用していた連結リストでは粒子の追加や消滅などを表現する際の扱いが容易という利点はあるが時間が経過するにつれて近傍の粒子がメモリアドレス上で散逸していき GPU で求められるコアレスアクセスに比べて効率が低くなるそこで連結リストから通常の配列に移行した前後の粒子を指し示すポインタを排除し粒子位置に応じて粒子を配列内でソートすることで粒子をコアレスにアクセスさせるようにしたソートのコストは生じるもののメモリアクセス最適化による高速化及び連結リストでのポインタの計算が不要になった効果の方が大きい GPUDirectRDMA[14] 異なるノード間にてGPU 同士で通信する場合 GPUDirect を用いることでネットワークのチップをGPU から直接経由させて途中 CPU にデータが達すること無く受信するGPU がネットワークのチップからデータを直接受け取る事ができるこの直接通信はソフトハード共に仕様の要求要件があるので詳細は [1,2] 等を参照して欲しい -39- 我々は kinaco に含まれる袖通信を全て GPUDirect にて実装することで CPU で通信した場合とほぼ同等の通信速度となる事を確認した 7. まとめ気象気候海洋モデルのGPU 対応について GPU 並びに気象気候海洋モデルの基礎を含めて今後の課題をまとめた気象気候海洋モデルのコード特性並びに開発者利用者の背景から今後の気象気候海洋モデルの開発は OpenACC とCUDA を併用した実装が最適と判断した将来の並列化案として気象システム内の要素を並列化させる事も検討していく OpenACC とCUDA の併用で気象気候海洋モデルに求められる条件 ( 高解像度化高速化精微化 ) を満たすにはCUDA の知識と経験が必須である我々はこれまで行ったCUDA による非静力学海洋モデルのGPU 実装をさらに推し進め気象気候海洋モデルを始めとする数値モデルのGPU 実装と高度化に役立つ知見を蓄積し次世代のエクサスケールスーパーコンピュータにて成果を出していくことを目指している参考文献 1.CUDACProgrammingGuide,NVIDIA 2.CUDACBestPracticeGuide,NVIDIA 3.John Cheng, Max Grossman, Ty McKercher,CUDAC プロフェッショナルプログラミング, インプレス 4.OpenACC,htps:/ 5.Satoh,M.,et al.(2014).the Nonhydrostatic Icosahedral Atmospheric Model:description and development. Progress in Earth and Planetary Science1(1):18. 6.Shimokawabe,T.,etal.(2010).An80 FoldSpeedup,15.0TFlopsFulGPU

12 Acceleration of Non-Hydrostatic Weather Model ASUCA Production Code. Proceedings of the 2010 ACM/IEEE InternationalConference for High Performance Computing, Networking, Storage and Analysis, IEEEComputerSociety: Yashiroetal.(2015),A Simulationof GlobalAtmosphereModelNICAM on TSUBAME2.5UsingOpenACC,GTC 2015,SanJose. 8.Consortium forsmal-scalemodeling, htp:/ 9.Fuhrer,O.,etal.(2018).Near-global climatesimulationat1km resolution: establishing a performance baseline on 4888GPUs with COSMO 5.0. Geosci.ModelDev.11(4): Posey S.and AdieJ.(2019),Sunny Skies Ahead! Versioning GPU acceleratedwrf to3.7.1,gtc 2019, SanJose. 11. Matsumura,Y.andHasumi,H (2008). A non-hydrostaticoceanmodelwitha scalable multigrid Poisson solver. OceanModeling24(1-2): Yamagishi, T. and Matsumura, Y. (2016),GPU Acceleration ofa Nonhydrostatic Ocean Model with a Multigrid Poisson/Helmholtz Solver, ProcediaComputerScience,80, Yamagishi, T. and Matsumura, Y. (2016).GPU Acceleration ofa Non- Hydrostatic Ocean Model with Lagrangian Particle Tracking. SupercomputingConference Yamagishi,T.,etal.(2018),AnMPI- CUDA Acceleration for a Non- Hydrostatic Ocean Model with GPUDirect RDMA Transfers, GTC Japan

NVIDIA Tesla K20/K20X GPU アクセラレータアプリケーションパフォーマンステクニカルブリーフ

NVIDIA Tesla K20/K20X GPU アクセラレータアプリケーションパフォーマンステクニカルブリーフ K20 GPU2 個に対するスピードアップ NVIDIA は Fermi アーキテクチャ GPU の発表によりパフォーマンスエネルギー効率の両面で飛躍的な性能向上を実現しハイパフォーマンスコンピューティング (HPC) の世界に変革をもたらしましたまた実際に GPU

RIST ニュース No.64(2018) GPU の気象 気候 海洋モデルへの適用とその最適化を目指した取り組み ActivitiestowardsanapplicationofGPU tonumerical weather/climate/oceanmodelanditsoptimization

RIST ニュース No.64(2018) GPU の気象気候海洋モデルへの適用とその最適化を目指した取り組み ActivitiestowardsanapplicationofGPU tonumerical weather/climate/oceanmodelanditsoptimization