Python + Numba CUDA で GPU プログラミングしてみるシリーズ
- その０： Python + Numba CUDA で GPU プログラミング入門以前
- その１： Python + Numba CUDA で画像処理 [← イマココ]
- その２： Python + Numba CUDA でライフゲーム
前回は GPU マシンに Numba と CUDA Toolkit を入れて、とりあえず GPU で簡単な計算をさせるプログラムを動かすことができた。
今回は CUDA のプログラミングモデルを理解しながら、GPU で画像を扱えるところまでやっていく。
CUDA のプログラミングモデル
- ゼロから始めるGPU Computing | G-DEP
このページがとても詳しくてわかりやすいのでここを読んでもらったほうが早い。
ホストとデバイス
- CPU 側のことをホスト、GPU 側のことをデバイスと呼ぶ
- デバイス側で実行される処理をカーネル関数を呼ぶ
- 並列処理をする流れは
0. カーネル関数をデバイスに転送する
0. データをデバイスのメモリに転送する
0. カーネル関数を起動する
スレッド・ブロック・グリッド
並列処理の並列数はスレッド・ブロック・グリッドという単位で管理される。
-
スレッドはカーネル関数が実行される最小単位
- ひとつのスレッドはひとつの CUDA コア内で実行される
-
ブロックはスレッドをまとめた単位
- 1ブロックあたり最大で 512 スレッドをまとめられる
- ３次元的な表現で管理でき、例えば「8x8x8 の３次元で合計 512 スレッド」のようにまとめられる
- ひとつのブロックの処理はひとつの SM (ストリーミングマルチプロセッサ) 内で実行される
-
グリッドはブロックをまとめた単位
- 1グリッドあたり最大で 65536x65536 ブロックをまとめられる
- ひとつのグリッドはひとつのデバイス (GPU) 内で実行される
カーネル関数を起動する際は「１ブロックあたりのスレッド数」「１グリッドあたりのブロック数」を指定する必要がある。
Numba CUDA の使い方
ざっくり解説するが、詳しくは公式ドキュメント見て欲しい。
Numba for CUDA GPUs — Numba documentation
カーネル関数の定義
-
@cuda.jitデコレータをつけて関数を定義するとそれがカーネル関数になる。 - カーネル関数内で
cudaオブジェクトから現在のスレッドに関する情報が得られる- これにより、現在のスレッドでは入力データのどの要素を対象に処理をすればいいかを計算することができる
例えば、入力が二次元配列 (行列) の場合は以下のようにして処理対象の座標 x, y を計算する。
from numba import cuda
import numpy as np
@cuda.jit
def increment_by_one(arr):
"""行列の要素に 1 を足すカーネル関数"""
tx = cuda.threadIdx.x
ty = cuda.threadIdx.y
bx = cuda.blockIdx.x
by = cuda.blockIdx.y
bw = cuda.blockDim.x
bh = cuda.blockDim.y
x = tx + bx * bw
y = ty + by * bh
arr[x, y] += 1
また、これは Numba CUDA の便利メソッドを利用すると以下のように書き換えることができる。
from numba import cuda
import numpy as np
@cuda.jit
def increment_by_one(arr):
"""行列の要素に 1 を足すカーネル関数"""
x, y = cuda.grid(2)
arr[x, y] += 1
カーネル関数の起動
定義したカーネル関数に対して、griddim と blockdim を指定した上で関数を実行する。
griddim = 1 # 1グリッドあたりのブロック数
blockdim = 3, 3 # 1ブロックあたりのスレッド数
# 入力する行列
arr = np.array([[0, 1, 2], [3, 4, 5], [6, 7, 8]])
increment_by_one[griddim, blockdim](arr)
print(arr) #=> [[1, 2, 3], [4, 5, 6], [7, 8, 9]]
画像を読み込んで NumPy 配列に変換する
画像処理を行うための事前準備として、今回利用する画像ファイルを読み込めるようにしておく。
Python では、PIL で読み込んだ画像オブジェクトを np.array に渡すことで画像を NumPy 配列に変換することができる。
import numpy as np
from PIL import Image
img = Image.open('lenna.gif')
arr = np.array(img)
NumPy 配列の画像は、PIL の関数を使うことで PIL 画像に戻すことができる。
Image.fromarray(img_inverted)
画像の色を反転させる
簡単な例として、画像の色を反転させるだけの処理を GPU で行ってみる。
from numba import cuda
import numpy as np
from PIL import Image
@cuda.jit
def invert_color(img_in, img_out):
"""画像の色を反転させるカーネル関数"""
x, y = cuda.grid(2)
if x < img_in.shape[0] and y < img_in.shape[1]:
img_out[x, y] = 0xFF - img_in[x, y]
# 画像を読み込んで NumPy 配列に変換する
img = np.array(Image.open('lenna.gif'))
# 変換結果を保存する NumPy 配列を作る
img_inverted = np.zeros(img.shape, dtype=np.uint8)
# データを GPU に転送してカーネル関数を実行する
invert_color[img.shape, 1](img, img_inverted)
img_inverted の配列を画像として表示してみると、反転できていることがわかる。
畳み込み演算
もう少し複雑な計算の例として、ラプラシアンフィルタでの画像のエッジ検出を試してみる。
ラプラシアンフィルタでは、以下のカーネル行列1を用いて畳み込み演算を行う。
\begin{bmatrix}
1 & 1 & 1 \\
1 & -8 & 1 \\
1 & 1 & 1
\end{bmatrix}
from numba import cuda
import numpy as np
from PIL import Image
@cuda.jit
def convolution(kernel, img_in, img_out):
"""畳み込み演算を行うカーネル関数"""
x, y = cuda.grid(2)
if x < img_in.shape[0] and y < img_in.shape[1]:
val = 0
for i in range(kernel.shape[0]):
if x + i == 0 or x + i == img_in.shape[0]:
continue
for j in range(kernel.shape[1]):
if y + j == 0 or y + j == img_in.shape[1]:
continue
val += img_in[x+i-1, y+j-1] * kernel[i, j]
img_out[x, y] = val if val > 0 else 0
# 画像を読み込んで NumPy 配列に変換する
img = np.array(Image.open('lenna.gif'))
# 変換結果を保存する NumPy 配列を作る
img_filtered = np.zeros(img.shape, dtype=np.uint8)
# 畳み込み演算に用いるカーネル行列を定義する
kernel = np.array([
[ 1, 1, 1],
[ 1, -8, 1],
[ 1, 1, 1],
])
# 畳み込み演算を実行する
convolution[img.shape, 1](kernel, img, img_filtered)
フィルタ後の画像を表示すると、ちゃんとエッジが検出できていることがわかる。
まとめ
GPU を用いて画像の畳み込み演算を行うことができた。
この程度の処理であれば CPU でやった方が速い (GPU はデータ転送のオーバーヘッドが大きい) のだが、より計算量の多い処理を行う際には GPU の恩恵が得られるに違いない。
ひとまず今回はここまで。
-
CUDA のカーネル関数と畳み込み演算のカーネル行列で「カーネル」という言葉が被っていて紛らわしいが、今回の話ではこれらは別物と考えてもらってよい ↩