ステンシル計算を対象とした 大規模GPUクラスタ向け 自動並列化フレームワーク 野村 達雄1 丸山 直也1 遠藤 敏夫1 松岡 聡1,2 1.東京工業大学 2.国立情報学研究所

背景 GPUクラスタの普及 GPUクラスタの難点 専門的な知識を持った一部の人の利用に留まる。もったいな い！ GPU一枚でTFlopsのピーク性能 高スループットのデバイスメモリ アプリケーションによっては大幅な性能向上 GPUクラスタの難点 GPU専用のプログラミング言語 CPUコード、GPUコード、ノード間並列のためのコード 複雑なメモリモデル。CPUメモリ、GPUメモリ 専門的な知識を持った一部の人の利用に留まる。もったいな い！

研究目的 開発者に特別な知識がなくてもGPUクラスタの性能 を享受できるようにしたい 開発者にGPUや並列性(CUDA,MPI)を意識させない CPUをシリアルに使うようにプログラミング 手動で最適化されたGPUコード並のパフォーマンス

提案 GPUクラスタ向けのDomain Specificな自動並列化フ レームワーク ドメインを限定することで抽象的な記述、高度な最適化が可能 ターゲットはステンシル計算 科学計算で頻繁に用いられており、実アプリが豊富 GPU化した場合の性能向上が著しい C言語を入力としたステンシル計算のコンパイラ、ライブラリ 開発者はGPUや並列化などを意識せずにステンシルを記述 フレームワークがCUDA,MPIを使った並列コードを生成

成果 フレームワークを一部実装 フレームワークを使ってGPUクラスタ上で動くステンシ ル計算を記述 C言語のステンシル記述をGPUクラスタ向けに自動並列化 フレームワークを使ってGPUクラスタ上で動くステンシ ル計算を記述 三次元拡散方程式 記述量はCPUシリアルで記述した場合と同程度 性能はCPUシリアルコードの35倍、手動 (GPU+MPI)で記述した ものの65%程度、 CPUをシリアルに使うのと同程度の労力で35倍速い

ステンシル計算の概要 隣接点計算 流体計算の計算カーネルとして よく表れる 有限差分法 空間を離散化し、各点の値をタイ ムステップでアップデート アップデートを決まった回数、ま たは値が収束するまで繰り返す メモリへのアクセスパターンが 単純であり、GPUで性能向上し やすい 二次元流体計算の例

ステンシルアプリケーションの例 T : 時間 X,Y : 空間 for (t = 0.0; t < T; t += dt) { for (y = 0; y < Y; y++) { for (x = 0; x < X; x++) { new_f(x,y)= e1*old(y,x) + e2*old(y,x-1) + e3*old(y,x+1) + e4*old(y-1,x) + e5*old(y+1,x); } 5点ステンシル

GPUクラスタ向けの実装 GPUクラスタ 袖領域 MPI GPU CPU CUDA 空間を分割して 各ノードに割り当てる

実装の複雑さ CPU,GPU,MPIコード 計算は簡潔、並列化のためのコー ドが多い 最適化のためのコード 袖領域の交換のためのコード GPUとGPUは直接通信できない GPU->CPU -> CPU->GPUの手順 最適化のためのコード ステンシルとは本質的には関係な い部分が複雑 Shimokawabe (sc2010)

提案フレームワーク 提案フレームワーク C言語によるアーキテクチャ独立、 GlobalViewの記述 CPU,GPU,MPIコードの生成 Code 2 Code 静的な解析による最適化 動的なコード生成、 オートチューニング チェックポインティング 最適な計算リソースの自動的な 選択 Parallerization Target Processors .c 提案フレームワーク NVidia GPUs AMD GPUs Nehalem Opteron … a.out MPI …

ステンシル計算の記述 ステンシル計算は主に三つの要素で構成される グリッド（空間）定義 DeclareGrid3D(real, float); グリッドの一点のアップデート方法を記述したカーネル関数 隣接点から値を読み、格子点の値を計算する float ret = grid_get(g,0,-1,0) + grid_get(g,0,1,0); return ret / 2.0; グリッドに対するカーネル関数の呼び出し グリッド上のすべての格子点についてカーネル関数が適用される grid_update(g, kernel);

グリッド上のすべての点に対してカーネルを呼び出す 7点ステンシルアプリケーションの記述例 DeclareGrid3D(real, float); float average(int x, int y, int z, grid3d_real g) { float ret = grid_get(g,0,0,0) + grid_get(g,-1,0,0) + grid_get(g,1,0,0) + grid_get(g,0,-1,0) + grid_get(g,0,1,0) + grid_get(g,0,0,-1) + grid_get(g,0,0,1); return ret / 7.0; } void computation(float *inbuff, float *outbuff) { grid3d_real g = grid3d_real_new(N,N,N); grid_copyin(g, inbuff); for (int t = 0; t < T; t += dt) grid_update(g, average); grid_copyout(g, outbuff); 定数位置離れた点の値 求められた点の値 グリッド上のすべての点に対してカーネルを呼び出す

グリッド上のすべての点に対してカーネルを呼び出す 7点ステンシルアプリケーションの記述例 DeclareGrid3D(real, float); float average(int x, int y, int z, grid3d_real g) { float ret = grid_get(g,0,0,0) + grid_get(g,-1,0,0) + grid_get(g,1,0,0) + grid_get(g,0,-1,0) + grid_get(g,0,1,0) + grid_get(g,0,0,-1) + grid_get(g,0,0,1); return ret / 7.0; } void computation(float *inbuff, float *outbuff) { grid3d_real g = grid3d_real_new(N,N,N); grid_copyin(g, inbuff); for (int t = 0; t < T; t += dt) grid_update(g, average); grid_copyout(g, outbuff); 定数位置離れた点の値 求められた点の値 // ステンシルカーネルから生成されるコードの擬似コード __device__ float average(int x, int y, int z, grid3d_real g) { float ret = g->buff[index(x,y,z)] + g->buff[index(x-1,y,z)] + g->buff[index(x+1,y,z)] + g->buff[index(x,y-1,z)] + g->buff[index(x,y+1,z)] + g->buff[index(x,y,z-1)] + g->buff[index(x,y,z+1)]; return ret / 7.0; } グリッド上のすべての点に対してカーネルを呼び出す

グリッド上のすべての点に対してカーネルを呼び出す 7点ステンシルアプリケーションの記述例 //カーネル呼び出しから生成される擬似コード __global__ void run_kernel(grid3d_real g) { 　　　int x = getx(), y = gety(); // GPUのスレッドIDから計算 for (int z = 0; z < g->dimz; z++) { g->back_buff[index(x,y,z)] = kernel(x, y, z, index); } void grid_update(grid3d_real g) { mpi_exchange_boundary(g); // MPIで袖領域の交換 run_kernel<<<blocks, threads>>>(g); // GPUカーネルの呼び出し DeclareGrid3D(real, float); float average(int x, int y, int z, grid3d_real g) { float ret = grid_get(g,0,0,0) + grid_get(g,-1,0,0) + grid_get(g,1,0,0) + grid_get(g,0,-1,0) + grid_get(g,0,1,0) + grid_get(g,0,0,-1) + grid_get(g,0,0,1); return ret / 7.0; } void computation(float *inbuff, float *outbuff) { grid3d_real g = grid3d_real_new(N,N,N); grid_copyin(g, inbuff); for (int t = 0; t < T; t += dt) grid_update(g, average); grid_copyout(g, outbuff); 定数位置離れた点の値 求められた点の値 グリッド上のすべての点に対してカーネルを呼び出す

記述の制約 カーネル関数の制約 格子点から相対位置で定数距離にある点から値を読み込みできる カーネル関数以外の部分の制約 袖領域のサイズを検出できる 読み込み限定のため、隣接点を誤って書き込むことがない カーネル関数以外の部分の制約 格子点への読み書きはできない データは分散しているため、ランダムなアクセスを許すとパフォーマ ンスを落とす可能性が高い 初期値の設定や計算後のデータはまとめてバッファからR/W

評価 フレームワークの一部を実装し、評価 目的 評価方法 簡潔な記述でGPUクラスタ向けのコード生成出来ることを確認 手動で記述した場合との記述量、性能の比較 評価方法 流体計算の一部である三次元拡散方程式（7点ステンシル）を 記述 手動実装とフレームワーク利用 計算精度 float

コードの記述量 空白行やコメント行を含めたおおよその行数 CPUシリアルの場合と同程度の記述量、内容もほぼ同等 生成されたコードと手動のコードの行数の違いは補助関数やタイプ 定義など 手動(CPU) 110行 手動(MPI-GPU) 250行 フレームワーク 生成コード 340行

性能評価用マシンの構成 CPU Model Intel Core i7 920 Clock 2.67 GHz Cores 4 physical cores GPU Tesla C2050 1.15 GHz Device Memory 3 GB Compute Capability 2.0 CUDA Runtime Version 3.0.14 CUDA Driver Version 195.36.15 Host Memory DDR3 12GB Network Infiniband DDR 20Gb/s OS CentOS 5.3 MPI OpenMPI 1.4.1 MAXノード数 10 ノード

1ノードでの性能(256x256x256) 35倍

考察 無駄なインデックス計算（約23GFlops） ループの違い（約7GFlops） これらを修正すれば手動のコードと同程度の性能になる。 grid_get(g, -1, 0, 0), grid_get(g, 1, 0, 0); 自動生成では各点について絶対的なアドレスを計算している、手動では 相対 もともと計算量の少ないカーネルではアドレス計算がボトルネックにな る ループの違い（約7GFlops） 手動: stencil_kernel() { for (z) …} z方向のループがステンシルカーネルに含まれるため、変数の再利用などが可能 自動: for (z) { stencil_kernel(); } z方向でループを回してステンシルカーネルを呼び出す これらを修正すれば手動のコードと同程度の性能になる。

Shared Memory vs. Global Memory

複数ノードの性能(256x256x512) 70% 65% 2ノードの時、自動生成では手動の65%程度の性能 10ノード時、70%程度の性能 スケールの仕方はほぼ同じであり、差は単一GPUのコード

関連研究 言語 対象 抽象度 ノード間並列 GPU 既存のプログラムとの親和性 HPF HPC全般 ○ × OpenMP C,Fortran - HMPP Chapel Pastha Haskel ステンシル Kamilらの研究 Fortran 我々の研究 C GPU版ASUCA

まとめ GPUクラスタを対象としたステンシル計算のコード生成 フレームワークを提案し、実装と評価を行った 拡散方程式を完結な表現で記述し、自動で並列化される ことを確認 CPUをシリアルに使った場合と同程度の記述量 1ノードを使った場合は手動/CPUの3500%,手動/GPUの 60%,10ノードで70%程度

今後の課題 生成されるコードの最適化、チューニング 宣言的な境界条件の記述 並列IOを行うためのインターフェイスの提供 Fortranの対応 インデックス計算の簡略化、ループ、分岐の最適化 スレッドブロック最適なサイズの自動探索 袖領域の交換による遅延を隠蔽するための最適化 宣言的な境界条件の記述 並列IOを行うためのインターフェイスの提供 Fortranの対応 実アプリでの性能、記述力の評価 地震動シミュレーション（東大古村先生） 気象予測モデル、ASUCA（気象庁） CPUとの協調、またはCPUクラスタへの対応