本文面向有一定編程基礎(chǔ)的開發(fā)者，系統(tǒng)介紹卷積運(yùn)算的數(shù)學(xué)原理以及如何用 CUDA 編寫高效的 GPU 并行卷積核函數(shù)。

一、卷積運(yùn)算是什么？

卷積的本質(zhì)是用一個小窗口（過濾器）在輸入數(shù)據(jù)上滑動，對每個局部區(qū)域做加權(quán)求和。

以圖像處理為例，對一個 的輸入矩陣，使用 的過濾器進(jìn)行卷積，計算過程如下

過濾器在輸入上逐步滑動，每滑動一步就輸出一個值，最終得到一張新的"特征圖"（feature map）。

二、CUDA 并行實現(xiàn) 2D 卷積

2.1 為什么用 GPU？

卷積中每個輸出像素的計算相互獨(dú)立，天然適合 GPU 的大規(guī)模并行架構(gòu)。GPU 可以讓數(shù)千個線程同時工作，每個線程負(fù)責(zé)一個輸出像素的全部計算。

2.2 線程映射策略

CUDA 采用二維線程網(wǎng)格對應(yīng)二維輸出圖像：

線程網(wǎng)格（Grid）
┌─────────────────────────────────────┐
│  Block(0,0)  Block(1,0)  Block(2,0) │
│  Block(0,1)  Block(1,1)  Block(2,1) │  ←── 每個 Block 內(nèi)有 N×N 個線程
│  Block(0,2)  Block(1,2)  Block(2,2) │
└─────────────────────────────────────┘
         ?  一一對應(yīng)  ?
       輸出圖像像素矩陣

每個線程的輸出像素坐標(biāo)由以下公式計算：

int outCol = blockIdx.x * blockDim.x + threadIdx.x;   // 列坐標(biāo)
int outRow = blockIdx.y * blockDim.y + threadIdx.y;   // 行坐標(biāo)

2.3 核函數(shù)完整實現(xiàn)

/**
 * 2D 卷積 CUDA 核函數(shù)
 * @param N      輸入圖像（行主序展開的一維數(shù)組）
 * @param F      卷積核（(2r+1)×(2r+1)，行主序展開）
 * @param P      輸出圖像
 * @param r      卷積核半徑（核大小 = 2r+1）
 * @param width  圖像寬度
 * @param height 圖像高度
 */
__global__ void convolution_2D_kernel(
    float *N,
    float *F,
    float *P,
    int r,
    int width,
    int height
) {
    // 1. 計算當(dāng)前線程負(fù)責(zé)的輸出像素坐標(biāo)
    int outCol = blockIdx.x * blockDim.x + threadIdx.x;
    int outRow = blockIdx.y * blockDim.y + threadIdx.y;

    // 超出圖像范圍的線程直接退出
    if (outCol >= width || outRow >= height) return;

    // 2. 初始化累加器
    float pValue = 0.0f;

    // 3. 遍歷卷積核的每個元素（大小為 (2r+1) × (2r+1)）
    int kernelSize = 2 * r + 1;
    for (int fRow = 0; fRow < kernelSize; fRow++) {
        for (int fCol = 0; fCol < kernelSize; fCol++) {
            // 計算對應(yīng)的輸入像素坐標(biāo)（以輸出像素為中心向外擴(kuò)展）
            int inRow = outRow - r + fRow;
            int inCol = outCol - r + fCol;

            // 4. 邊界檢查：超出圖像范圍的位置視為 0（zero-padding 語義）
            if (inRow >= 0 && inRow < height && inCol >= 0 && inCol < width) {
                pValue += F[fRow * kernelSize + fCol]   // 卷積核元素（正確的線性索引）
                        * N[inRow * width + inCol];     // 輸入像素
            }
        }
    }

    // 5. 寫入輸出圖像
    P[outRow * width + outCol] = pValue;
}

注意：卷積核以 float * 傳入，必須用 F[fRow * kernelSize + fCol] 進(jìn)行線性尋址，不能寫 F[fRow][fCol]（后者在 float * 上無法編譯）。

2.4 Host 端調(diào)用

void launch_convolution(
    float *h_N, float *h_F, float *h_P,
    int r, int width, int height
) {
    int kernelSize = 2 * r + 1;
    size_t imgBytes    = width * height * sizeof(float);
    size_t kernelBytes = kernelSize * kernelSize * sizeof(float);

    // 分配 GPU 內(nèi)存
    float *d_N, *d_F, *d_P;
    cudaMalloc(&d_N, imgBytes);
    cudaMalloc(&d_F, kernelBytes);
    cudaMalloc(&d_P, imgBytes);

    // 拷貝數(shù)據(jù)到 GPU
    cudaMemcpy(d_N, h_N, imgBytes,    cudaMemcpyHostToDevice);
    cudaMemcpy(d_F, h_F, kernelBytes, cudaMemcpyHostToDevice);

    // 配置線程塊與網(wǎng)格（推薦 16×16，共 256 線程/Block）
    dim3 blockDim(16, 16);
    dim3 gridDim(
        (width  + blockDim.x - 1) / blockDim.x,
        (height + blockDim.y - 1) / blockDim.y
    );

    // 啟動核函數(shù)
    convolution_2D_kernel<<<gridDim, blockDim>>>(d_N, d_F, d_P, r, width, height);

    // 拷回結(jié)果
    cudaMemcpy(h_P, d_P, imgBytes, cudaMemcpyDeviceToHost);

    cudaFree(d_N);
    cudaFree(d_F);
    cudaFree(d_P);
}

三、控制發(fā)散

3.1 什么是控制發(fā)散？

GPU 以 Warp（32 個線程）為單位調(diào)度執(zhí)行。同一 Warp 內(nèi)的線程必須執(zhí)行相同的指令分支，否則不滿足條件的線程會被掛起等待，稱為控制發(fā)散（Control Divergence）。

在卷積中，邊緣像素的輸入窗口超出圖像范圍，觸發(fā) if (inRow >= 0 && ...) 分支的跳過，而中心線程正常執(zhí)行，導(dǎo)致同一 Warp 內(nèi)分支不一致：

同一 Warp 內(nèi)（32 個線程）
┌──────────────────────────────────────┐
│ 線程 0-3：邊緣像素，部分跳過 if 分支  │  ← 閑置等待
│ 線程 4-31：中心像素，正常執(zhí)行全部計算 │  ← 正常工作
└──────────────────────────────────────┘
            ↓ 實際效果
   Warp 效率下降（部分線程空轉(zhuǎn)）

3.2 實際影響有多大？

控制發(fā)散只出現(xiàn)在圖像邊緣，其嚴(yán)重程度與以下因素成反比：

以常見場景為例： 圖像 + 卷積核，邊緣像素約占總量的 0.3%，控制發(fā)散的影響可忽略不計。

結(jié)論：對于大圖像 + 小卷積核（深度學(xué)習(xí)典型場景），無需針對控制發(fā)散專門優(yōu)化。只有在處理小圖像或極大卷積核時，才值得考慮針對性方案。

四、進(jìn)一步優(yōu)化方向

本文介紹的是基礎(chǔ)實現(xiàn)。實際生產(chǎn)環(huán)境中，還有以下優(yōu)化手段：

4.1 共享內(nèi)存 Tiling

全局內(nèi)存（Global Memory）延遲高達(dá)數(shù)百個時鐘周期。相鄰線程會重復(fù)讀取輸入圖像的相同區(qū)域，可以將輸入 Tile 預(yù)加載到共享內(nèi)存（Shared Memory，延遲約為全局內(nèi)存的 1/30），大幅減少全局內(nèi)存訪問次數(shù)。

全局內(nèi)存訪問模式（基礎(chǔ)版）    共享內(nèi)存 Tiling（優(yōu)化版）
每個線程獨(dú)立讀取 (2r+1)2    每個 Block 協(xié)作加載一塊 Tile
個全局內(nèi)存位置               到共享內(nèi)存，再從共享內(nèi)存讀取
→ 大量重復(fù)讀取               → 重復(fù)讀取命中緩存，帶寬壓力大降

4.2 常量內(nèi)存緩存卷積核

卷積核 在整個計算過程中只讀不寫，且所有線程訪問完全相同的數(shù)據(jù)，非常適合存入常量內(nèi)存（Constant Memory）：

__constant__ float d_F[MAX_KERNEL_SIZE];   // 存入常量內(nèi)存

常量內(nèi)存有專用緩存，廣播讀取效率極高。

4.3 線程粗化（Thread Coarsening）

讓單個線程負(fù)責(zé)多個輸出像素，減少線程調(diào)度開銷和 Warp 同步消耗，在某些架構(gòu)上可提升吞吐量。

參考延伸閱讀：NVIDIA CUDA Programming Guide、Deep Learning by Ian Goodfellow（第 9 章：卷積網(wǎng)絡(luò)）