Kaggle Notebook で始める CUDAプログラミング入門
はじめに
こんにちは！kaerururu です。
普段、私は PyTorch や cuML といったライブラリを通じて、実装の詳細を意識することなくCUDA の恩恵を受けています。CUDA プログラミングができるようになるともっと速いコードが書けそうです。
この記事では、NVIDIA公式のドキュメントなどを頼りにCUDAプログラミングを学んだ記録を共有します。C++での実装から、PythonのラッパーライブラリであるNumbaを使った実装まで、具体的なコードと実行速度の比較まで実施しました。
この記事で作成したコードは、Kaggle Notebook で公開しています。
CUDAとは？
CUDA (Compute Unified Device Architecture) は、NVIDIAが開発したGPU向けの並列コンピューティングプラットフォーム兼プログラミングモデルです。
CUDAプログラミングでは、主にC++言語を拡張してコードを記述します。CPUで実行される部分（ホストコード）と、GPUで並列実行される関数（カーネル関数）を明確に分けて記述するのが特徴です。
- ホスト (Host): CPUとそのメモリのこと。プログラム全体の制御や、GPUへの処理の依頼を行います。
- デバイス (Device): GPUとそのメモリのこと。ホストから依頼された計算処理を、大規模な並列実行で高速に処理します。
-
カーネル (Kernel):
__global__というキーワードで宣言される、デバイス（GPU）上で実行されるC++関数のことです。このカーネルが、多数のスレッドによって同時に実行されます。
計算の流れは、一般的に以下のようになります。
- ホスト（CPU）からデバイス（GPU）へデータを転送する。
- ホストがカーネル関数を呼び出し、デバイスに計算処理を命令する。
- デバイスが何千ものスレッドを使って並列に計算を実行する。
- ホストがデバイスから計算結果のデータを受け取る。
実行環境
今回は、誰でも無料でGPUが使えるGoogle ColabとKaggle Notebook (T4×2) の2つの環境で実験を行いました。それぞれの環境と主要なライブラリのバージョンは次の通りです。
| 環境 | Python |
nvidia-smi (CUDA Ver) |
nvcc -V |
g++ --version |
GPU |
|---|---|---|---|---|---|
| Kaggle Notebook | 3.11.13 | 12.6 | 12.5.82 | 11.4.0 | NVIDIA T4 |
| Google Colab | 3.12.12 | 12.4 | 12.5.82 | 11.4.0 | NVIDIA T4 |
※バージョンは実行時期によって変動する可能性があります。
実験
今回は、シンプルな例を題材に CUDA プログラミングに入門します。同じ長さの2つの配列の各要素を足し合わせる処理です。
y[i] = x[i] + y[i]
この単純な処理を、(1) 素のC++ (CPU)、(2) CUDA C++ (GPU)、(3) Numba (Python on GPU) の3つの方法で実装し、その速度を比較します。
1. C++ (CPU) での実装【ベースライン】
まずは、比較対象となるCPUでの実装です。forループを使って、配列の要素を一つずつ順番に計算します。
#include <iostream>
#include <math.h>
#include <chrono>
// function to add the elements of two arrays
void add(int n, float *x, float *y)
{
for(int i = 0; i < n; i++)
y[i] = x[i] + y[i];
}
int main(void)
{
int N = 1<<20; // 1M elements
float *x = new float[N];
float *y = new float[N];
// initialize x and y arrays on the host
for(int i=0; i<N; i++) {
x[i] = 1.0f;
y[i] = 2.0f;
}
// Run on 1M elements on the CPU
auto start = std::chrono::high_resolution_clock::now();
add(N, x, y);
auto end = std::chrono::high_resolution_clock::now();
std::chrono::duration<double, std::milli> elapsed = end - start;
std::cout << "CPU execution time: " << elapsed.count() << " ms\n";
// Check for errors (all values should be 3.0f)
float maxError = 0.0f;
for(int i=0; i<N; i++)
maxError = fmax(maxError, fabs(y[i]-3.0f));
std::cout << "Max error: " << maxError << std::endl;
// Free memory
delete [] x;
delete [] y;
return 0;
}
コンパイルして実行します。-O3オプションでコンパイラの最適化を有効にしています。
g++ -O3 add.cpp -o add
./add
2. CUDA C++ (GPU) での実装
GPU上で並列実行するカーネル関数 add_gpu を定義します。
スレッド、ブロック、グリッド
CUDAでは、計算を実行する最小単位を スレッド (Thread) と呼びます。このスレッドが、カーネル関数のコードを実行します。そして、スレッドは ブロック (Block) というグループにまとめられ、さらにブロックは グリッド (Grid) という単位で管理されます。
-
スレッド: 1つの計算を実行する最小単位。
threadIdx.xでブロック内の自身のIDを取得できます。 -
ブロック: スレッドの集まり。同じブロック内のスレッドは、高速な共有メモリを使って連携できます。
blockIdx.xでグリッド内の自身のIDを、blockDim.xでブロック内のスレッド数を取得できます。 - グリッド: ブロックの集まり。カーネルを呼び出す際の最小単位です。
これらの階層構造を使い、blockIdx, blockDim, threadIdx という組み込み変数から、各スレッドが担当すべき配列のインデックスを一意に計算します。
int i = blockIdx.x * blockDim.x + threadIdx.x;
この式がCUDAプログラミングの基本中の基本です。
#include <iostream>
#include <math.h>
#include <cuda_runtime.h> // For CUDA Event
// Kernel function to add the elements of two arrays
__global__
void add_gpu(int n, float *x, float *y)
{
int index = blockIdx.x * blockDim.x + threadIdx.x;
int stride = blockDim.x * gridDim.x;
for (int i = index; i < n; i += stride)
y[i] = x[i] + y[i];
}
int main(void)
{
int N = 1<<20; // 1M elements
// --- GPU Timer Setup ---
cudaEvent_t start, stop;
cudaEventCreate(&start);
cudaEventCreate(&stop);
// Allocate Unified Memory -- accessible from CPU or GPU
float *x, *y;
cudaMallocManaged(&x, N*sizeof(float));
cudaMallocManaged(&y, N*sizeof(float));
// initialize x and y arrays on the host
for(int i=0; i<N; i++) {
x[i] = 1.0f;
y[i] = 2.0f;
}
// Prefetch the x and y arrays to the GPU
cudaMemPrefetchAsync(x, N*sizeof(float), 0, 0);
cudaMemPrefetchAsync(y, N*sizeof(float), 0, 0);
// Recommended parallel pattern on T4 GPU
int blockSize = 256;
int numBlocks = (N + blockSize - 1) / blockSize;
cudaEventRecord(start);
add_gpu<<<numBlocks, blockSize>>>(N, x, y);
cudaEventRecord(stop);
cudaEventSynchronize(stop);
// Calculate GPU time
float milliseconds = 0;
cudaEventElapsedTime(&milliseconds, start, stop);
// Check for errors (all values should be 3.0f)
float maxError = 0.0f;
for(int i=0; i<N; i++)
maxError = fmax(maxError, fabs(y[i]-3.0f));
std::cout << "Max error: " << maxError << std::endl;
std::cout << "GPU Kernel time: " << milliseconds << " ms" << std::endl;
// Free memory
cudaFree(x);
cudaFree(y);
// Destroy events
cudaEventDestroy(start);
cudaEventDestroy(stop);
return 0;
}
カーネル呼び出し構文
kernel_name<<<Grid, Block>>>(arguments);
<<<...>>> の部分で、グリッド数とブロックあたりのスレッド数を指定します。今回は、1ブロックあたり256スレッドとし、全要素をカバーするのに必要なブロック数を計算してグリッド数としています。これはGPUの性能を引き出すための一般的な設定の一つです。
コンパイル方法
CUDAコードは、NVIDIAが提供する nvcc (NVIDIA CUDA Compiler) を使ってコンパイルします。
nvcc add_gpu.cu -o add_gpu -arch=sm_75
./add_gpu
オプションとして -arch=sm_75 を指定しました。これは「どのGPUアーキテクチャ向けにコードをコンパイルするか」を指定するものです。KaggleやColabで使われている NVIDIA T4 の場合、Compute Capabilityが 7.5 (Turingアーキテクチャ) なので、sm_75 です。デフォルト値が sm_75 なので T4 GPU を使用する場合、指定の必要はないですが、GoogleColab ではこの指定がないと、GPUが正しく認識されずコードが動かなかったです。
GPUのCompute Capabilityは、こちらのNVIDIAのページで確認できます。
deviceQuery
グリッド数とブロックあたりのスレッド数の最大値はどのように知ることができるのでしょうか。cuda-samples リポジトリにある deviceQuery を実行するとシステム内に存在する CUDA デバイスのプロパティを列挙することができます。リポジトリのソースコードをコンパイルして実行ファイルを作成します。今回私の実験コードでは、%%writefile で必要なソースコードのみをカレントディレクトリ内に配置しています。(本文ではコードは省略)
%%writefile helper_string.h
...
%%writefile helper_cuda.h
...
%%writefile deviceQuery.cpp
...
%%writefile CMakeLists.txt
...
!SMS=75 cmake .
!make
コンパイルできたら実行します。
!./deviceQuery
実行結果は次のとおりです。
Kaggle の T4×2 Notebook は文字通り T4 が 2つ接続されているので、Device 0 と Device 1 に "Tesla T4" が認識されています。
./deviceQuery 実行結果
./deviceQuery Starting...
CUDA Device Query (Runtime API) version (CUDART static linking)
Detected 2 CUDA Capable device(s)
Device 0: "Tesla T4"
CUDA Driver Version / Runtime Version 12.6 / 12.5
CUDA Capability Major/Minor version number: 7.5
Total amount of global memory: 15095 MBytes (15828320256 bytes)
(040) Multiprocessors, (064) CUDA Cores/MP: 2560 CUDA Cores
GPU Max Clock rate: 1590 MHz (1.59 GHz)
Memory Clock rate: 5001 Mhz
Memory Bus Width: 256-bit
L2 Cache Size: 4194304 bytes
Maximum Texture Dimension Size (x,y,z) 1D=(131072), 2D=(131072, 65536), 3D=(16384, 16384, 16384)
Maximum Layered 1D Texture Size, (num) layers 1D=(32768), 2048 layers
Maximum Layered 2D Texture Size, (num) layers 2D=(32768, 32768), 2048 layers
Total amount of constant memory: 65536 bytes
Total amount of shared memory per block: 49152 bytes
Total shared memory per multiprocessor: 65536 bytes
Total number of registers available per block: 65536
Warp size: 32
Maximum number of threads per multiprocessor: 1024
Maximum number of threads per block: 1024
Max dimension size of a thread block (x,y,z): (1024, 1024, 64)
Max dimension size of a grid size (x,y,z): (2147483647, 65535, 65535)
Maximum memory pitch: 2147483647 bytes
Texture alignment: 512 bytes
Concurrent copy and kernel execution: Yes with 3 copy engine(s)
Run time limit on kernels: No
Integrated GPU sharing Host Memory: No
Support host page-locked memory mapping: Yes
Alignment requirement for Surfaces: Yes
Device has ECC support: Enabled
Device supports Unified Addressing (UVA): Yes
Device supports Managed Memory: Yes
Device supports Compute Preemption: Yes
Supports Cooperative Kernel Launch: Yes
Supports MultiDevice Co-op Kernel Launch: Yes
Device PCI Domain ID / Bus ID / location ID: 0 / 0 / 4
Compute Mode:
< Default (multiple host threads can use ::cudaSetDevice() with device simultaneously) >
Device 1: "Tesla T4"
CUDA Driver Version / Runtime Version 12.6 / 12.5
CUDA Capability Major/Minor version number: 7.5
Total amount of global memory: 15095 MBytes (15828320256 bytes)
(040) Multiprocessors, (064) CUDA Cores/MP: 2560 CUDA Cores
GPU Max Clock rate: 1590 MHz (1.59 GHz)
Memory Clock rate: 5001 Mhz
Memory Bus Width: 256-bit
L2 Cache Size: 4194304 bytes
Maximum Texture Dimension Size (x,y,z) 1D=(131072), 2D=(131072, 65536), 3D=(16384, 16384, 16384)
Maximum Layered 1D Texture Size, (num) layers 1D=(32768), 2048 layers
Maximum Layered 2D Texture Size, (num) layers 2D=(32768, 32768), 2048 layers
Total amount of constant memory: 65536 bytes
Total amount of shared memory per block: 49152 bytes
Total shared memory per multiprocessor: 65536 bytes
Total number of registers available per block: 65536
Warp size: 32
Maximum number of threads per multiprocessor: 1024
Maximum number of threads per block: 1024
Max dimension size of a thread block (x,y,z): (1024, 1024, 64)
Max dimension size of a grid size (x,y,z): (2147483647, 65535, 65535)
Maximum memory pitch: 2147483647 bytes
Texture alignment: 512 bytes
Concurrent copy and kernel execution: Yes with 3 copy engine(s)
Run time limit on kernels: No
Integrated GPU sharing Host Memory: No
Support host page-locked memory mapping: Yes
Alignment requirement for Surfaces: Yes
Device has ECC support: Enabled
Device supports Unified Addressing (UVA): Yes
Device supports Managed Memory: Yes
Device supports Compute Preemption: Yes
Supports Cooperative Kernel Launch: Yes
Supports MultiDevice Co-op Kernel Launch: Yes
Device PCI Domain ID / Bus ID / location ID: 0 / 0 / 5
Compute Mode:
< Default (multiple host threads can use ::cudaSetDevice() with device simultaneously) >
> Peer access from Tesla T4 (GPU0) -> Tesla T4 (GPU1) : Yes
> Peer access from Tesla T4 (GPU1) -> Tesla T4 (GPU0) : Yes
deviceQuery, CUDA Driver = CUDART, CUDA Driver Version = 12.6, CUDA Runtime Version = 12.5, NumDevs = 2
Result = PASS
この辺りが、グリッド数とブロックあたりのスレッド数の上限に関係しています。
grid size は非常に大きい値が記載されていますが、これらは論理アドレス空間として用意してあるだけで、全てが利用されることはないことに注意です。
Maximum number of threads per multiprocessor: 1024
Maximum number of threads per block: 1024
Max dimension size of a thread block (x,y,z): (1024, 1024, 64)
Max dimension size of a grid size (x,y,z): (2147483647, 65535, 65535)
3. Numba (Python) での実装
C++のコードに不慣れでも Numba を使用すれば、Python で CUDA プログラミングができます。
Numbaは、PythonコードをJIT (Just-In-Time) コンパイルすることで高速化するライブラリで、CUDAもサポートしています。@cuda.jit というデコレータを付けるだけで、Python関数をCUDAカーネルに変換できます。
import numpy as np
from numba import cuda
@cuda.jit
def add_numba(x, y):
i = cuda.grid(1)
if i < x.size:
y[i] = x[i] + y[i]
def main():
n = 1 << 20
x = np.ones(n, dtype=np.float32)
y = np.full(n, 2.0, dtype=np.float32)
d_x = cuda.to_device(x)
d_y = cuda.to_device(y)
threads_per_block = 256
blocks_per_grid = (n + (threads_per_block - 1)) // threads_per_block
start_event = cuda.event()
end_event = cuda.event()
start_event.record()
add_numba[blocks_per_grid, threads_per_block](d_x, d_y)
end_event.record()
end_event.synchronize()
elapsed_time_ms = cuda.event_elapsed_time(start_event, end_event)
print(f"Numba GPU execution time: {elapsed_time_ms} ms")
y = d_y.copy_to_host()
absolute_errors = np.abs(y - 3.0)
max_error = np.max(absolute_errors)
print(f"{max_error=}")
if __name__ == "__main__":
main()
C++版と比較して、メモリ管理（cudaMallocManaged など）がNumpyライクな cuda.to_device, copy_to_host で直感的に書けます。カーネルの呼び出し方も似ています。
パフォーマンス比較と考察
Kaggle Notebook (NVIDIA T4) で、約100万要素の配列の加算を実行した結果は以下のようになりました。
| 実装方法 | 実行時間 (ms) | CPUに対する高速化倍率 |
|---|---|---|
| C++ 最適化 (CPU) | 0.356 ms |
約1.0x |
| CUDA C++ (GPU) | 0.054 ms |
約6.6x |
| Numba (GPU) | 0.048 ms |
約7.4x |
※時間は実行のたびに多少変動します。
考察
今回の実験では次のようになりました。
- CUDA 実装の実行時間はコンパイラの最適化オプションを指定した CPU 実装の 6~7 倍
- 配列計算の部分のみの実行時間のみでみると、CUDA C++ と Numba の実行時間はほぼ同速
CUDAカーネルのみの実行時間は C++ も Numba もほぼ同速 (誤差の範囲) ですが、全体の実行時間で見るとメモリコピーなどのオーバヘッド分 C++ の方が速かったです。
しかし、C++のコードを書く手間やコンパイルの手間を考えると、Pythonのコードに数行追加するだけでこの速度が得られるNumba は非常に使いやすそうです。
まとめ
今回は、CUDAプログラミングの第一歩として、配列の加算を題材にC++とPython (Numba) での実装とパフォーマンス比較を行いました。
- CUDAはC++を拡張し、GPUの並列処理能力を最大限に引き出すための強力なフレームワークであること。
- カーネル関数、スレッド、ブロック、グリッドといった基本概念。
-
nvccを使ったコンパイルと、-archオプションの紹介。 - Numba を使えば、Pythonから手軽にCUDAのパワーを利用できること。
これらを手を動かしながら学ぶことができました。今回は簡単な題材でしたが、次回はもっと実践的な題材で CUDA プログラミングの理解を深めたいです。
Discussion