GPU Convolutional Layer：CUDA Kernel 與 GEMM 化

重點總覽 (Overview)

卷積層的可平行性 (Parallelism in the Convolutional Layer)

回顧 16-Deep-Learning/02-Convolutional-Neural-Networks-Layers 的 minibatch forward 程式碼 (Fig. 16.12) 共有 7 層迴圈：

for n          // minibatch sample      ┐
  for m        // output feature map    │  4 層 "easy" parallel
    for h      // output row            │  iterations = N·M·H_out·W_out
      for w    // output col            ┘
        for c      // input channel     ┐  3 層需 atomic 累加
          for p    // filter row        │  (read-modify-write 同一個 Y)
            for q  // filter col        ┘  → 預設保持序列

• 公式：可用平行 iteration 總數 = N · M · H_out · W_out，其中 H_out = H - K + 1、W_out = W - K + 1。
• 這四層已足夠餵滿 GPU，故 c/p/q 留在 thread 內序列執行，避免 atomic 的開銷。

直接 CUDA 推論 Kernel (Direct CUDA Inference Kernel)

Thread → data 映射設計

• 每個 thread 算 1 個 output pixel。
• 2D thread block：每個 block 負責一個 output feature map 上 TILE_WIDTH × TILE_WIDTH 的 tile（TILE_WIDTH=16 → 256 threads/block）。
• 3D grid (M, T, N)：
  • X 維 = M（output feature map 索引）。
  • Z 維 = N（minibatch sample 索引）。
  • Y 維 = T = H_grid · W_grid：把垂直/水平 tile 索引線性化塞進單一維度（因 X、Z 已被佔用）。

grid (M, T, N)，block (TILE_WIDTH, TILE_WIDTH, 1)

blockIdx.x = m  ──► 哪一張 output feature map
blockIdx.z = n  ──► 哪一個 minibatch sample
blockIdx.y      ──► row-major 線性化的 tile 編號

  一張 output feature map (W_grid=2, H_grid=2):
   ┌─────────┬─────────┐
   │ tile0,0 │ tile0,1 │   blockIdx.y = 0,1,2,3
   │ (by=0)  │ (by=1)  │   (row-major)
   ├─────────┼─────────┤
   │ tile1,0 │ tile1,1 │   tile 內由 threadIdx.x/y
   │ (by=2)  │ (by=3)  │   定位 16×16 pixel
   └─────────┴─────────┘

還原 tile 座標：
  vertical_tile   = blockIdx.y / W_grid
  horizontal_tile = blockIdx.y % W_grid

Host code (Fig. 16.13)

#define TILE_WIDTH 16
W_grid = W_out / TILE_WIDTH;     // horizontal tiles per output map
H_grid = H_out / TILE_WIDTH;     // vertical   tiles per output map
int T  = H_grid * W_grid;        // tiles per map, linearized into grid.y
dim3 blockDim(TILE_WIDTH, TILE_WIDTH, 1);
dim3 gridDim(M, T, N);
ConvLayerForward_Kernel<<<gridDim, blockDim>>>(/* C, W_grid, K, X, W, Y */);

Kernel (Fig. 16.15)

__global__ void
ConvLayerForward_Kernel(int C, int W_grid, int K,
                        float* X, float* W, float* Y) {
    int m = blockIdx.x;                                        // output feature map
    int h = (blockIdx.y / W_grid) * TILE_WIDTH + threadIdx.y;  // output row
    int w = (blockIdx.y % W_grid) * TILE_WIDTH + threadIdx.x;  // output col
    int n = blockIdx.z;                                        // minibatch sample
    float acc = 0.0f;
    for (int c = 0; c < C; c++)        // sum over all input channels
        for (int p = 0; p < K; p++)    // KxK filter
            for (int q = 0; q < K; q++)
                acc += X[n, c, h + p, w + q] * W[m, c, p, q];  // row-major linearize!
    Y[n, m, h, w] = acc;
}

效能定位

GEMM 化：im2col Unrolling (Formulating Convolution as GEMM)

中心思想（Chellapilla et al., 2006）：展開並複製 input pixel，讓計算某一 output pixel 所需的所有輸入剛好排成 unrolled 矩陣的一整欄 (column)，於是整個 forward 變成一次大矩陣乘法，丟給 cuBLAS GEMM。

三個矩陣的形狀

• 高度 C·K·K：每個 output element 來自 C 張 input map、各取 K·K 個元素。
• 寬度 H_out·W_out：output feature map 的像素數。
• 注意 M 不影響 X_unroll——所有 output map 共用同一份展開輸入。

範例 (Fig. 16.16): C=3, 每張 3×3；M=2，每張 output 2×2；filter 2×2

  filter matrix        ×     X_unroll           =    Y
  (2 × 12)                  (12 × 8...實為 4)        (2 × 4)
 ┌────────────┐         ┌──────────────┐        ┌──────────┐
 │ map0 weights│        │col=各 output │        │ out map0 │
 │ map1 weights│        │  pixel 的    │   →    │ out map1 │
 └────────────┘         │  輸入 patch  │        └──────────┘
                        └──────────────┘
 一欄展開示意 (output (0,0))：
   X 三張 map 的 patch 串接 → (1,2,1,1, 0,2,0,3, 1,2,0,1)ᵀ
   filter 串接             → (1,1,2,2, 1,1,1,1, 0,1,1,0)
   內積 = 14 = Y[0,0,0]

Expansion ratio（膨脹比）

由於 convolution 的 patch 彼此重疊，input pixel 被重複多次：

3×3 input map 中：
  4 個角 pixel  → 各用 1 次  (4·1)
  4 個邊 pixel  → 各用 2 次  (2·4)
  1 個中心pixel → 用 4 次    (4·1)
  展開後 = 4 + 8 + 4 = 16 個，原本 9 個 → ratio = 16/9 ≈ 1.8

• 一般公式：
  • expansion_ratio = (C·K·K·H_out·W_out) / (C·H_in·W_in)
  • 當 feature map 遠大於 filter 時 → 趨近 K·K。
• 範例：(3·2·2·2·2)/(3·3·3) = 16/9。

Unroll CUDA kernel (Fig. 16.18)

設計：每個 thread 從一張 input map 蒐集 K·K 個元素，填入 unrolled 矩陣某欄的一段；總 thread 數 = C · H_out · W_out，用 1D block。

__global__ void
unroll_Kernel(int C, int H, int W, int K, float* X, float* X_unroll) {
    int t = blockIdx.x * blockDim.x + threadIdx.x;
    int H_out = H - K + 1, W_out = W - K + 1;
    int W_unroll = H_out * W_out;            // width of unrolled matrix
    if (t < C * W_unroll) {
        int c       = t / W_unroll;          // which input channel
        int w_unroll= t % W_unroll;          // column = linearized output pixel
        int h_out   = w_unroll / W_out;
        int w_out   = w_unroll % W_out;
        int w_base  = c * K * K;             // section start row for channel c
        for (int p = 0; p < K; p++)
            for (int q = 0; q < K; q++) {
                int h_unroll = w_base + p*K + q;          // row in unrolled matrix
                X_unroll[h_unroll, w_unroll] = X[c, h_out + p, w_out + q];
            }
    }
}

執行流程與記憶體取捨

• 資料佈局：input N×C×H×W、filter 預先排成 M×C×K²、output M×(H_out·W_out)。
• 記憶體足跡控制：不一次展開整個 minibatch，只配置 1 個 X_unroll[C·K·K, H_out·W_out] buffer，逐 sample 重用 → 每個 sample 先 unroll_Kernel、再 GEMM。

考試/面試重點 (Exam / Test Patterns)

Related Notes

• 16-Deep-Learning/02-Convolutional-Neural-Networks-Layers
• 16-Deep-Learning/04-cuDNN-Library-and-Summary
• 07-Convolution/03-Tiled-Convolution-and-Halo-Handling
• 03-Multidimensional-Grids-And-Data/04-Matrix-Multiplication-Kernel
• 05-Memory-Architecture-And-Data-Locality/02-Tiling-and-Tiled-Matrix-Multiplication
• 06-Performance-Considerations/01-Memory-Coalescing