Mathematica
による
GPU
計算入門
--
手持ちの
MacBookPro(2-core 4-thread)
を使い、
GPU
プログラミングを
(
実演
)
。
＠Mathematica研究会第4回 2012.3.10 松田裕幸
はじめに
Mathematica8.0からGPUプログラミング/実行環境 CUDA(Common Unified Device Architecture)とOpenCLをサポートす
るCUDALinkが提供された。 CUDALinkは従来のMathLink方式ではなく、MathematicaカーネルはCUDA用にコンパイ
ル・ライブラリ化されたDLL(Dynamic Link Libray)を直接呼び出す方式を採用。
本来であれば、この実装によりMathematicaで高速なGPU計算環境が実現できるはずなのだが。。。。 Mathematicaの
ドキュメントにはGPUプログラミングの高速化技術が例によっていくつか紹介されているが、実際には、そうした主
としてメモリ管理に関する高速化より、CUDALink自身のプロトコル管理のオーバヘッドが大きすぎ、Mathematicaと
の間の通信がほとんどなく、GPU上で綺麗に並列化され、かつ長時間走る問題しか高速化は望めない。
ただし。GPUプログラミングがどんなものであるか、あるいはGPUプログラミングを経験するのにMathematicaは超
お手頃。理由は３つ。
1) カーネルコードに集中できる。
2) インタラクティブプログラミングーデバッグ、開発が容易 3) 教育効果
MathematicaからCUDALinkが使えるまで
GPUに対するCUDA用ドライバのインストール
特定のGPU用専用ドライバを間違えていストールしないように。汎用版を利用。
http://developer.nvidia.com/cuda-toolkit-40
Mathematica用CUDA Resourceパッケージのインストール
Version 8.0.4.1 BuildNumber -1
QualifiedName CUDAResources-OSX-8.0.4.1
Description {ToolkitVersion -> 4.0, MinimumDriver -> 270.0} MathematicaVersion 8.0.1+
SystemIDs MacOSX-x86-64
Location /Users/yuko/Library/Mathematica/Paclets/Repository/CUDAResources-OSX-8.0.4.1
PC内蔵GPUエンジン
GPUエンジンスペック CUDAプロセッサーコア 48
ギガフロップス 182
オープンGL 2.1
In[7]:=
Out[7]//TableForm=
Name GeForce GT 330M Maximum Grid Dimensions 65535
65535 1 Maximum Block Dimensions 512
512 64 Maximum Threads Per Block 512 Maximum Shared Memory Per Block 16384 Total Constant Memory 65536 Multiprocessor Count Total Memory 3219259392 Warp Size 32 Core Count 48
GPUメモリ管理
GPUプログラミングのイメージ
Mathematicaコード
1000*1000スレッドを作成、開始、指標を各スレッドに渡し、結果のInputDataを関数funで計算し、結果をOutputData
に置く。
スレッドという表現はGPUアーキテクチャから来る発想で、プログラム的にはカーネルという表現が使われる。 CUDAカーネルプログラムのファイル拡張子は *.cu を使用。
アーキテクチャ
CUDAプログラムにおけるスレッドは、GPU上のGrid/Block/Threadにマップされる。二次元の画像処理の場合だと、
一次元セルラオートマトンの場合
上図はCUDAプログラムの一般的なサイクルの詳細を示している．
Allocate GPU Memory -- GPUのメモリを割り当てる．GPUメモリはCPUメモリとは別のもので，プログラマは割当
てコピーを管理しなければならない．
Configure Threads -- スレッドの構成を設定する．問題に適したブロックとグリッドの次元を選ぶ．
Launch Threads -- 設定したスレッドを開始する
CUDA Kernel Code
Synchronize Threads -- CUDAスレッドを同期させると，デバイスがGPUメモリ上でその後の操作を行う前にすべて
のタスクを確実に完了させるようにできる．
Copy GPU Memory to CPU -- スレッドが完了したら，メモリはGPUからCPUにコピーして返される．
Free GPU -- GPUメモリが解放される．
CUDALinkは”CUDA Kernel Code”以外の面倒をすべてみてくれる。逆にいえば、CUDALinkのアーキテクチャの出来
如何によってCUDAコードの実行性能は大きく変わってくる。そして残念ながら現状(Mathematica 8.0.4)でもあまり改
善は見られない。
メモリの種類
大域メモリ Global memory
GPU上で利用できる，最も豊富だが最も遅いメモリ．このメモリは128，256，512MBのパッケージングがある．すべ
てのスレッドは大域メモリの要素にアクセスできるが，パフォーマンス上の理由からこれらのアクセスは最小限に抑 えられ，さらに制約がある．
テクスチャメモリ Constant Texture
テクスチャメモリは大域メモリと同じ場所にあるが，読取り専用である．テクスチャメモリは大域メモリに見られる パフォーマンスの低下はないが，サポートしている型はchar，int，floatのみである．
コンスタントメモリ Constant Memory
グリッド内のどのスレッドからもアクセスできる高速なコンスタントメモリ．メモリはキャッシュされるが，全部で
64KBに制限される．
共有メモリ Shared Memory
特定のブロックに局所的な高速メモリ．現在のハードウェア上では共有メモリ量はブロックごとに16KBに制限され
ている．
局所メモリ Local Memory (in thread = in kernel)
使用するGPUの計算能力を確認する
CUDALink/tutorial/Programming
In[8]:=
In[9]:=
Out[9]=
計算能力は，デバイスが行える操作を表す．現段階では計算能力は1.1，1.2，1.3，2.0のみがある．主な違いを以下に
リストする．
計算能力 計算能力 1.0 基本実装 1.1 アトミック演算 1.2 共有アトミック演算とWarp Vote関数 1.3 倍精度演算のサポート 2.0 倍精度，L2キャッシュ，並行カーネル
CUDAカーネルを書く
__device__ float f(float x) { return tanf(x);
}
__global__ void secondKernel(float * diff, float h, mint listSize) { mint index = threadIdx.x + blockIdx.x * blockDim.x;
float f_n = f(((float) index) / h); float f_n1 = f((index + 1.0) / h); if( index < listSize) {
diff[index] = (f_n1 - f_n) / h; }
}
__global__
下線2個 GPU上で実行される（C関数として実行）ことを指示する関数修飾子。
__device__
__global__関数、あるいは他の__device__関数から呼び出される関数の関数修飾子。C関数から直接呼び出されること
はない。
void
CUDAプログラミングでは __global__関数は値を返さず、出力は関数引数で指定したバッファへの書き込みによって
行う。この例では、arryが入力と出力を兼ねている。
型 (float * diff, float h, mint listSize)
参考：CUDALink/ref/CUDAFunctionLoad 詳細部分に記載。左端がMathematica側でコンパイルする時に用いる関数プ
ロトタイプの型。真ん中がCUDAカーネルで使用する型。なおmintはMathematica固有の型。
_Integer mint Mathematica integer Integer32 int 32-bit integer _Real Real_t GPU real type Double double machine double Float float machine float
Local | Shared mint local or shared memory parameter {Local | Shared, type} mint local or shared memory parameter
mint index = threadIdx.x + blockIdx.x * blockDim.x;
ドキュメントには以下のように書かれているがわかりにくい。
threadIdx — 現在のスレッドの指標．スレッド指標は0とblockDim - 1の間である
blockIdx — 現在のブロックの指標．ブロック指標は0とgridDim - 1の間である
blockDim — ブロックサイズ次元．
gridDim — グリッドサイズ次元．
リストサイズlistSizeを256、コンパイル時指定するブロックサイズを16とすると、 threadIdx.x = 0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15, 0, 1,2, ....
blockIdx.x = 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1.. blockDim.x = 16
となる。
float f_n = f(((float) index) / h);
関数fの呼び出し。変数index, h共に整数なので除算による切り捨てを防ぐためindexの型をfloatに変換。
if( index < listSize)
出力用バッファ領域を超えて書き込みが起きないようガード。
実行例
CUDAカーネルコード
実際にCUDAカーネルコードを書くときはエディタ等でファイルに保存し、それをMathematicaからImortして使うの
が普通。
In[13]:=
Out[14]=
ここで{16,16}は16*16二次元ブロック数を表す。上記 blockDimがこれに相当。
またコンパイル時、”Defines”によってカーネルコード中で使用した環境変数に対し定数を与えることができる。
In[15]:=
実行
In[21]:=
Out[24]=
CUDA アプリケーション
CUDALink/tutorial/Applications
Image Binarize
カーネルコード
In[25]:=
コンパイル
In[26]:=
入力、出力バッファの準備。入力バッファはCPU側、出力バッファはCUDAMemoryAllocateによってGPU側に確保。
In[27]:=
2値化スレッシュフォールドを150として実行。
In[31]:=
Out[31]=
結果が出力バッファoutputに入っているのでそれを表示。
In[32]:=
Out[32]=
Mathematica組み込み関数で同じ結果を確認。
In[33]:=
GPU側に確保したメモリを開放。
In[34]:=
CUDA関数
CUDALink/tutorial/Functions
MathematicaにはCUDA化された関数がいくつか用意されている。現状ではあまり使い物にならないと思うが。。。
ヒストグラム
In[35]:=
ヒストグラム用データの収集。計算を速くするためにMathematicaコードをコンパイルしている。
In[36]:=
Out[37]=
表示。
In[38]:=
CUDAFunctionLoad
CUDALink/ref/CUDAFunctionLoad
CUDAプログラミングを行う際に必要な細かい話、および少し大きめのサンプルコードが多数紹介されている。
一次元セルラオートマトン
ソースコード。Mathematicaをインストールし、CUDALinkパッケージをダウンロードするとCUDALinkフォルダの直
下にサンプルコードが格納される。
In[39]:=
コンパイル。GPUにバイナリロード。
In[41]:=
バッファを用意。
ルール30を128回実行。
In[46]:=
結果を表示。
In[47]:=
Out[47]=
