Convolution 練習題 (Practice - Convolution Fundamentals and the Basic Parallel Kernel)

Convolution：每個輸出元素是「對應輸入元素 + 中心對稱鄰居」的 weighted sum，權重由 filter 給定。大小 2r+1 為奇數，使加權和對稱於被計算的中心元素，每側剛好 r 個鄰居（r = radius）。書中為避免與 CUDA kernel function 撞名，把權重陣列一律稱 filter (F)。

Ghost cells 是靠近陣列端點時、落在邊界外不存在的輸入元素；多數應用補預設值 0（如 audio 假設錄音前後音量為 0；也有應用填最近的有效邊緣值）。2D 更複雜，因為輸出可能同時碰到水平、垂直或兩者皆有的邊界（例如 P_{1,0} 缺 2 行 + 1 列）。

用 2D grid，每 thread 算一個 output pixel：outCol = blockIdx.x*blockDim.x + threadIdx.x、outRow = blockIdx.y*blockDim.y + threadIdx.y（與 Ch.3 ColorToGrayscale 相同）。瓶頸是 global memory bandwidth：inner loop 2 ops / 載入 8 bytes（N 4B + F 4B）= 0.25 OP/B。次要問題是邊緣 thread 處理不同數量 ghost cells 造成的 control divergence（大圖小 filter 時影響輕微到可忽略）。

宣告：__constant__ float F[..][..];，必須放在所有函式之外（global scope）。Host 拷貝用專用 API cudaMemcpyToSymbol(dest, src, size)（不是 cudaMemcpy），它告知 runtime 這份資料 kernel 期間不會被改。容量約 64 KB。Kernel 把 F 當全域變數直接存取，不必當參數傳入（跨檔需 extern 宣告）。

三特性：(1) 體積小（radius 通常 ≤ 7，3D filter ≤ 7³=343）；(2) kernel 執行中內容不變（read-only）；(3) 所有 thread 以相同順序存取 F。當一個 warp 內 32 個 thread 同時讀同一位址（F 索引與 threadIdx 無關），硬體可一次廣播，提供龐大頻寬。因 F 視為永遠命中、不耗 DRAM，OP/B 由 0.25 翻倍到 0.5 OP/B（2 ops / 4 bytes，只剩 N 的 4B）。

公式：IN_TILE_DIM = OUT_TILE_DIM + 2*FILTER_RADIUS（每方向各延伸 r 以涵蓋 halo）。Halo = 鄰居 tile 的真實內部元素；ghost = 陣列外、補 0 的不存在元素。方案 A 第一個 if 管 ghost cell（邊界外補 0），第二個 if 管 deactivation（block 比 output tile 大，停用外圈 FILTER_RADIUS 層 threads）。__syncthreads() 放在載入完整 input tile 之後、開始計算之前（barrier）。

用 y[i] = F[0]·x[i-1] + F[1]·x[i] + F[2]·x[i+1]，邊界 ghost = 0：

• y[0] = 2·0 + 1·4 + 4·1 = 8
• y[1] = 2·4 + 1·1 + 4·3 = 21
• y[2] = 2·1 + 1·3 + 4·2 = 13
• y[3] = 2·3 + 1·2 + 4·3 = 20
• y[4] = 2·2 + 1·3 + 4·0 = 7

輸出 = {8, 21, 13, 20, 7}。

令 M = 2r+1、IN_TILE_DIM = T + 2r = T + M - 1。

• (a) blocks = ⌈N/T⌉²（每 block 一個 output tile）。
• (b) threads/block = IN_TILE_DIM² = (T+M-1)²（方案 A：block 尺寸 = input tile）。
• (c) shared memory = (T+M-1)² × 4 bytes。
• (d) 5×5 (r=2)、IN=32 → OUT=28：ratio = 28²·25·2 / (32²·4) = 39200/4096 ≈ 9.57 OP/B（上界 (2r+1)²·2/4 = 12.5）。

兩者把複雜度在「載入 ↔ 計算」之間搬家。方案 A：載入單純（每 thread 載 1 元素）、但需停用外圈 threads且 input≠output tile（不易為 2 的次方、divergence 較多）；halo 已在 shared，halo DRAM 流量為 0。方案 B：input = output = TILE_DIM（可為 2 的次方、divergence 較少、載入只需檢查一般邊界），但計算迴圈較複雜（內部讀 shared、halo 讀 global 且仍要判 ghost），halo 倚賴 L2 cache 命中（機率高但非保證）。

(a) 透明度不同：cache（含 constant cache）對程式透明，程式只存取原始 global 變數、硬體自動快取並記住位址；shared memory 是 scratchpad，須手動宣告 __shared__ 並明確拷貝。(b) filter 越大，每個 N 元素被約 (2r+1)² 個 thread 重用，故上界 (2r+1)²·2/4 越高；但 halo 隨 r 增大、把 output tile 擠小，使實際 ratio 與上界落差變大。小 tile 危險：8×8 input tile 配 5×5 filter 僅 3.13 OP/B（halo 佔比過高），3D convolution 因 on-chip memory 不足常被迫用小 tile，reduction 效果大打折扣。