Lecture 7: GPU 架构与 CUDA 编程（GPU Architecture & CUDA Programming）（日期：Oct 14, 2025）

目录 · ← l6 · l8 →

Lecture 7: GPU 架构与 CUDA 编程（GPU Architecture & CUDA Programming）（日期：Oct 14, 2025）

概述：本讲从 GPU 的历史讲起——GPU 原本是为实时 3D 游戏渲染设计的专用处理器，后来人们发现它对”大规模数据上执行相同计算”（数据并行）极其擅长，于是 2007 年 NVIDIA 随 Tesla 架构推出了 CUDA，让 GPU 以通用计算模式运行任意程序。本讲的核心目标是：① 掌握 CUDA 编程抽象（thread / block / grid 两级线程层级、host/device 分离、shared/global memory、__syncthreads、原子操作）；② 理解这些抽象在现代 GPU（以 NVIDIA V100 为例）上是如何实现的（SIMT、warp、SM、线程块调度器）。课程反复强调的思考题是：CUDA 究竟是数据并行模型、共享地址空间模型还是消息传递模型？它与 ISPC 的 gang/task、pthreads 有何异同？

注意：本讲对应 Assignment 3: A Circle Renderer in CUDA（CUDA 圆形渲染器）——你需要用本讲的全部 CUDA 概念（block/grid 配置、shared memory、__syncthreads、scan 等）在 GPU 上实现一个高性能圆形渲染器。

一、核心概念与定义

1. CUDA（Compute Unified Device Architecture）

定义：NVIDIA 于 2007 年随 Tesla 架构推出的”C 风格”编程语言与运行时，用于在 GPU 的 compute mode（通用计算模式）硬件接口上编写程序。它相对底层：CUDA 的抽象与当代 GPU 的能力/性能特征非常贴近，设计目标是保持低抽象距离（low abstraction distance）。
现实类比：CUDA 之于 GPU，就像 C 之于 CPU——不给你铺满玫瑰花的抽象（比如自动并行化），而是把硬件的真实结构（线程层级、共享内存、屏障）直接暴露给你，让你自己安排，换来得天独厚的性能控制力。
公式/图示：无；CUDA 是语言+运行时。

2. Kernel（内核，`global` 函数）

定义：用 __global__ 修饰、在 GPU（device）上以 SPMD 方式执行的函数。一次 kernel launch（kernel<<<grid, block>>>(args)）会批量启动（bulk launch）成千上万个 CUDA 线程，每个线程执行同一份 kernel 代码，但通过内置变量 threadIdx / blockIdx / blockDim 区分自己处理的数据。
现实类比：kernel 就像一张”施工图纸”（函数体），而批量启动就是复印几万份图纸发给几万个工人，每个工人按自己胸牌上的编号（threadIdx/blockIdx）负责不同工位。注意每个 worker 的程序计数器是独立的——这是与 ISPC gang 的本质区别之一（ISPC 是编译期把实例编译成 SIMD 指令，CUDA 是运行时动态判断）。

公式/图示：

一次 kernel 启动：  myKernel<<<numBlocks, threadsPerBlock>>>(args);
                   └─ 启动 numBlocks × threadsPerBlock 个 CUDA 线程

3. CUDA 线程（CUDA thread）与线程层级：thread / block / grid

定义：CUDA 线程是逻辑控制流（与 pthread 的抽象类似，但实现天差地别——见后文）。线程按两级层级组织：若干线程组成一个 thread block（线程块），若干 block 组成一个 grid（网格）。线程 ID 最高可以是 3 维的（dim3），方便处理天然 N 维的问题。
现实类比：想象一家工厂：grid 是整个车间（一批订单），block 是班组（一个班组内的工人必须同时在场、能互相递工具——对应共享内存与 __syncthreads），thread 是工人。订单（block）之间互相独立，车间调度员想先做哪单就先做哪单。

公式/图示：

Grid（本次 kernel 启动的全部线程，由 numBlocks 指定）
├── Block (0,0) ── 12 个线程：T(0,0) T(1,0) T(2,0) T(3,0)
│                              T(0,1) T(1,1) ...（2D 情形）
├── Block (1,0) ── 12 个线程
└── Block (2,0) ── 12 个线程
线程全局坐标：  i = blockIdx.x * blockDim.x + threadIdx.x
                j = blockIdx.y * blockDim.y + threadIdx.y

例：dim3 threadsPerBlock(4,3); dim3 numBlocks(3,2); → 6 个 block × 12 线程 = 72 个 CUDA 线程。

4. Host / Device 与分布式地址空间（distributed address space）

定义：CUDA 程序被静态地分成两半：host 代码（普通 C/C++，串行跑在 CPU 上）与 device 代码（kernel，跑在 GPU 上）。host 与 device 拥有不同的地址空间：host memory 与 device global memory 之间靠 cudaMalloc / cudaMemcpy（如 cudaMemcpyHostToDevice）搬运数据。在 host 端直接解引用 device 指针是非法的。
现实类比：就像两座隔海的城市，不能直接隔海递包裹，必须走”港口+轮船”（cudaMemcpy）。这正是课程前几讲讲的分布式内存（消息传递）风格的地址空间——cudaMemcpy 让你联想到什么？它跟 MPI 的 memcpy 一样在”两个地址空间之间”移动数据。

公式/图示：

Host（CPU）                      Device（GPU）
┌─────────────────┐    cudaMemcpy    ┌──────────────────────┐
│ Host memory     │◄────────────────►│ Device global memory │
│ 地址空间         │                  │ （DRAM，所有线程可读写）│
└─────────────────┘                  └──────────────────────┘

5. Device 内存模型：private / shared / global memory

定义：kernel 内部可见三种地址空间：per-thread private memory（每线程私有，通常是寄存器）、per-block shared memory（块内共享，片上高速存储）、per-program global memory（所有线程可见，位于 DRAM）。三种地址空间反映了程序中不同粒度的局部性——这是 GPU 高效实现 CUDA 的关键（如果预先知道某些线程访问同一批变量，调度器就能把共享数据放到高速片上存储）。
现实类比：global memory 是”公共仓库”（慢但容量大），shared memory 是”本班组共用的工作台”（快，但只限本班组），private memory 是”个人口袋里的笔记本”（最快）。

公式/图示：见第 4 点的图，补充三层结构：

Device 内部：
┌──────────────────────────────────────────┐
│ Global memory（所有 block、所有线程共享）    │
│   ┌────────────────────────────────────┐  │
│   │ Shared memory（block 内所有线程共享） │  │
│   │   ┌──────────────────────────────┐ │  │
│   │   │ Private memory（线程私有，寄存器）│ │  │
│   │   └──────────────────────────────┘ │  │
│   └────────────────────────────────────┘  │
└──────────────────────────────────────────┘

6. Warp（线程束）与 SIMT（Single Instruction, Multiple Thread）

定义：warp 是 32 个连续的 CUDA 线程组成的一组（block 内线程 0-31 为 warp 0，32-63 为 warp 1……）。GPU 硬件把同一个 warp 内 32 个线程的指令流取出后，动态检查它们是否执行同一条指令；若是，就用 SIMD ALU 让 32 个线程同时执行——NVIDIA 称之为 SIMT。warp 不属于 CUDA 编程模型，而是现代 NVIDIA GPU 上重要的实现细节。若 warp 内线程指令不一致（如 if 分支两侧都执行），则发生分支发散（divergent execution），性能受损。
现实类比：warp 像一列 32 节车厢的火车，共用一个火车头（取指/译码单元）；只要所有车厢目的地一致（同一条指令），火车头就能一次把整列车拖过去。哪节车厢想”变道”（走不同分支），整列车就得分成两趟跑。

公式/图示：

一个 256 线程的 block → 8 个 warp（256 / 32）
V100 一个 sub-core：最多可调度/交织 16 个 warp
warp 指令执行：16 个 fp32 ALU 跑 32 线程的指令 → 需要 2 个时钟

7. Shared Memory（共享内存，`shared`）

定义：kernel 内用 __shared__ 声明的、按 block 分配的片上高速存储（V100 上 shared + L1 合计 128 KB/SM）。块内所有线程可读写；块与块之间不共享。用途：把 global memory 中会被多次复用的数据一次性搬进片上，避免反复访问 DRAM。
现实类比：做菜时把冰箱（DRAM）里要用的所有食材一次性搬到厨房台面（shared memory）上，后面每个步骤都在台面上取料，而不是每用一次就跑一趟冰箱。
公式/图示：__shared__ float support[THREADS_PER_BLK + 2]; —— 每个 block 一份，块内线程共享。

8. `__syncthreads()`（块内屏障）与原子操作（atomic）

定义：__syncthreads() 是 block 内所有线程必须到达的屏障（barrier）：谁先到谁等待，全部到齐才继续。CUDA 还提供原子操作（如 atomicAdd(float* addr, float amount)），可作用于 global 和 shared 地址；另有 host/device 同步（kernel 返回时隐含对所有线程的屏障）。
现实类比：接力赛的”交接区”——所有队员必须都到达交接区，下一棒才能一起出发。少一个队员（比如某些线程没调用 __syncthreads()），全队卡死。

公式/图示：

线程 0: 写 shared ─┐
线程 1: 写 shared ─┼─► __syncthreads() ─► 所有线程才允许读 shared
线程 2: 写 shared ─┘

9. SM（Streaming Multiprocessor）与 sub-core

定义：SM 是 GPU 上的一个”多线程 SIMD 核”，V100 芯片上有 80 个 SM；每个 SM 由 4 个 sub-core 组成，每个 sub-core 有自己的 warp selector、取指/译码单元和一组 SIMD 功能单元（16 个 fp32 mul-add、16 个 int、8 个 fp64、load/store 单元、tensor core 单元）。寄存器堆共 256 KB/SM（每 sub-core 64 KB），分给最多 64 个 warp。
现实类比：SM 是一家”车间”，sub-core 是车间里的 4 条”装配线”，每条装配线同一时刻只伺候一个 warp，但可以在 16 个 warp 之间快速切换（交织执行）来隐藏延迟。

公式/图示：V100 关键数字：

80 SMs × 4 sub-cores × 16 fp32 ALUs = 5,120 fp32 mul-add ALUs = 12.7 TFLOPs
（mul-add 计为 2 flops）
最多 80 × 64 = 5,120 个并发 warp = 163,840 个并发 CUDA 线程/芯片
L2 cache 6 MB；HBM 16 GB，带宽 900 GB/sec（4096-bit 接口）

10. 线程块调度器（Thread Block Scheduler）与工作分配

定义：CUDA 的核心假设是：thread block 之间没有依赖，可以按任意顺序执行。GPU 上的硬件工作调度器（work scheduler）把 block（”工作”）按动态调度策略映射到 SM 上，只要满足资源约束（每 block 需要的线程上下文数、shared memory 字节数）。block 一旦完成，其资源（shared memory、warp 上下文）立即释放给下一个 block。注意 CUDA 程序里没有 num_cores 这个概念——同一个 kernel 可以不加修改地跑在 6 核或 16 核的 GPU 上（类似数据并行模型里的 forall）。
现实类比：公司接了一批订单（block），调度台按各产线的空闲情况随时派单，先完成先释放产线；派单员不关心订单之间的先后（因为订单之间本来就独立）。

公式/图示：

Grid（1000 个 block）→ GPU Work Scheduler → Core 0 / Core 1（fictitious 双核 GPU）
Step 1: host 发送 kernel 启动命令（EXECUTE convolve, NUM_BLOCKS=1000）
Step 2: 调度器把 block 0 映射到 core 0（预留 128 线程上下文 + 520B shared）
Step 3: 继续映射 block 1、2、3……（交错映射）
       核心容量：每 core 只能驻留 2 个 block（3 × 520B > 1.5KB shared）
Step 4: block 0 完成 → 资源释放
Step 5: block 4 调度到 core 0 ……（依次类推，共 1000 个）

11. 内存带宽（Memory Bandwidth）

定义：GPU 从 DRAM 读写数据的速率。CPU 时代内存带宽约 150-300 GB/sec（DDR5），高端 GPU 约 900 GB/sec-1 TB/sec（HBM，4096-bit 接口）。带宽是 GPU 程序最重要的性能天花板之一：shared memory 之所以快，是因为它是片上存储，不走 DRAM 带宽。
现实类比：带宽像”港口吞吐量”——每秒能从仓库（DRAM）运多少货到车间。计算再快，货（数据）运不进来也是白搭（这就是后续讲的 bandwidth-bound）。
公式/图示：带宽 = 总线宽度 × 时钟频率 × 每时钟传输次数（如 4096-bit × 1.4 GHz × 2 ≈ 900+ GB/s）。

12. GPGPU（General-Purpose computation on GPU）与 CUDA 的诞生

定义：2002-2003 年研究者用”hack”方式把 GPU 当作数据并行机器用：把输出图像大小设成数组大小（如 512×512），画两个恰好盖满屏幕的三角形，让 fragment shader 对每个像素执行一次计算——shader 函数被 map 到 512×512 个元素上。2004 年斯坦福图形学实验室的 Brook 语言把 GPU 抽象成流处理器（stream + kernel）。2007 年 NVIDIA Tesla 架构提供首个非图形专用的 compute mode 接口：分配 buffer、上传 kernel 二进制、launch(myKernel, N) 以 SPMD 方式运行 N 个实例——这比图形接口 drawPrimitives() 简单得多，CUDA 由此诞生。
现实类比：GPU 本来是”只会画画的专用打印机”，GPGPU 时代人们发现”只要把计算画成图像”就能让它算任何东西；CUDA 则干脆给打印机装上了”通用打印”按钮。

公式/图示：

hack 方法： 图像 512×512 ←→ 数组 512×512
            fragment shader（纯函数，跑在每个像素上） ←→ 对每个数组元素执行 f

二、代码示例与详细解说

示例 1：`matrixAdd`——2D block/grid 配置与 `<<<>>>` 启动

// matrixAdd.cu —— 编译：nvcc -o matrixAdd matrixAdd.cu（需 NVIDIA GPU 与 CUDA Toolkit）
#include <cuda_runtime.h>
#include <cstdio>

const int Nx = 12;
const int Ny = 6;

// ============ device 代码：kernel 定义（运行在 GPU 上） ============
__global__ void matrixAdd(float A[Ny][Nx], float B[Ny][Nx], float C[Ny][Nx])
{
    int i = blockIdx.x * blockDim.x + threadIdx.x;   // 列坐标（全局）
    int j = blockIdx.y * blockDim.y + threadIdx.y;   // 行坐标（全局）
    C[j][i] = A[j][i] + B[j][i];
}

// ============ host 代码：串行执行在 CPU 上 ============
int main()
{
    float *A, *B, *C;   // host 指针
    // （完整程序还需 cudaMalloc 设备内存并用 cudaMemcpy 拷贝数据，见示例 2 说明）
    dim3 threadsPerBlock(4, 3);                     // 每个 block 4×3 = 12 个线程
    dim3 numBlocks(Nx / threadsPerBlock.x,          // 网格 3×2 = 6 个 block
                   Ny / threadsPerBlock.y);
    // 这次启动共创建 72 个 CUDA 线程：6 个 block，每个 12 个线程
    matrixAdd<<<numBlocks, threadsPerBlock>>>(A, B, C);
    cudaDeviceSynchronize();                        // host/device 同步
    return 0;
}

【代码做了什么？】

host 端用 dim3 声明 threadsPerBlock(4,3) 与 numBlocks(3,2)（Nx/threadsPerBlock.x = 12/4 = 3，Ny/threadsPerBlock.y = 6/3 = 2）。
matrixAdd<<<numBlocks, threadsPerBlock>>>(A,B,C) 是一次 bulk launch：一次性创建 72 个 CUDA 线程（6 block × 12 线程），每个线程独立执行 kernel 函数体。
kernel 内每个线程用公式 i = blockIdx.x * blockDim.x + threadIdx.x、j = blockIdx.y * blockDim.y + threadIdx.y 计算自己在 12×6 矩阵中的全局坐标，然后执行 C[j][i] = A[j][i] + B[j][i]。
该调用是同步语义：host 端调用返回时所有线程都已终止（隐含全线程屏障）。

【并行机制解说】

线程如何创建：不是像 pthread_create 那样逐个创建（那要分配栈、OS 控制块），而是一次启动一整网格的线程，成本极低。
工作如何分配：程序里没有 num_cores——block 是”工作单元”，由 GPU 硬件调度器动态映射到任意数量的 SM 上，因此同样的程序能跑在 6 核或 16 核 GPU 上（对应概念：thread block scheduler、数据并行 forall 精神）。
数据如何共享：此例无共享，三个数组都在 global memory。
同步点：kernel 返回时的隐式屏障；本 kernel 块间无依赖，调度顺序任意。
对应概念：thread/block/grid 层级、kernel 启动语法、host/device 分离。

示例 2：1D 卷积 v1——每输出元素一个线程（朴素版）

// convolve_v1.cu —— 编译：nvcc -o convolve_v1 convolve_v1.cu
#include <cuda_runtime.h>

#define THREADS_PER_BLK 128

// kernel：对长度为 N 的输入做 3 点滑动平均：output[i] = (input[i]+input[i+1]+input[i+2]) / 3
__global__ void convolve(int N, float* input, float* output)
{
    int index = blockIdx.x * blockDim.x + threadIdx.x;  // 线程局部变量：全局线程号
    float result = 0.0f;                                // 线程局部变量：累加器
    for (int i = 0; i < 3; i++)
        result += input[index + i];                     // 3 次 global memory 读
    output[index] = result / 3.f;                       // 1 次 global memory 写
}

// host 代码
int main()
{
    int N = 1024 * 1024;                    // 约 100 万个输出元素
    float *devInput, *devOutput;
    cudaMalloc(&devInput,  sizeof(float) * (N + 2));  // 输入数组（多 2 个边界元素）
    cudaMalloc(&devOutput, sizeof(float) * N);        // 输出数组
    // （此处应使用 cudaMemcpy(devInput, hostInput, ..., cudaMemcpyHostToDevice)
    //   初始化 devInput 的内容，代码从略）
    // 启动 N/128 = 8192 个 block，共 1M 个线程，每个线程算 1 个输出元素
    convolve<<<N / THREADS_PER_BLK, THREADS_PER_BLK>>>(N, devInput, devOutput);
    cudaDeviceSynchronize();
    return 0;
}

【代码做了什么？】

每个 CUDA 线程负责一个输出元素：先由 index 定位自己是第几个输出，然后串行执行 3 次加法累加 input[index..index+2]，最后写回 output[index]。
这是典型的 “one thread per output element”（每个输出元素一个线程）的数据并行分解：1M 个输出 → 1M 个线程 → 8192 个 block。
注意 input 在 host 端无法直接访问：devInput 是设备地址空间指针，host 只能通过 cudaMemcpy 与之交互（对应概念：分布式地址空间）。

【并行机制解说】

工作分配：每个线程独立处理一个输出，线程之间零通信、零同步——这是最容易并行化的模式，GPU 上数万并发线程轻而易举。
效率问题（本讲的伏笔）：每个线程需要 3 次 global load（input[index]、[index+1]、[index+2]），且相邻线程的访问窗口互相重叠——相邻线程会重复读取同一个输入元素。整个 kernel 共执行 3×128 次 load per block，而输入数据总量只有 130 个元素/block。这是巨大的带宽浪费，引出示例 3。
对应概念：memory bandwidth 意识、数据并行分解（one thread per element）。

示例 3：1D 卷积 v2——用 shared memory 暂存输入（本讲重点）

// convolve_v2.cu —— 编译：nvcc -o convolve_v2 convolve_v2.cu
#include <cuda_runtime.h>

#define THREADS_PER_BLK 128

// kernel：与 v1 相同功能，但先把 block 需要的输入区间搬进 shared memory
__global__ void convolve(int N, float* input, float* output)
{
    __shared__ float support[THREADS_PER_BLK + 2];  // 每 block 分配：128 个 + 2 个 halo
    int index = blockIdx.x * blockDim.x + threadIdx.x;  // 线程局部变量

    // 1) 所有线程协作把本 block 的"支撑区"从 global 载入 shared memory
    support[threadIdx.x] = input[index];
    if (threadIdx.x < 2) {                            // 前 2 个线程额外加载 2 个 halo 元素
        support[THREADS_PER_BLK + threadIdx.x] =
            input[index + THREADS_PER_BLK];
    }

    __syncthreads();                                  // 2) 块内屏障：全部载入完成后才继续

    float result = 0.0f;                              // 线程局部变量
    for (int i = 0; i < 3; i++)
        result += support[threadIdx.x + i];           // 3) 从 shared memory 读，不再访问 global

    output[index] = result / 3.f;                     // 4) 写结果到 global memory
}

// host 代码与 v1 完全相同：
//   cudaMalloc 分配 devInput(N+2)、devOutput(N)；
//   convolve<<<N/THREADS_PER_BLK, THREADS_PER_BLK>>>(N, devInput, devOutput);

【代码做了什么？】

步骤 1：块内 128 个线程协作加载：每个线程把 input[index] 写入 support[threadIdx.x]；前 2 个线程（threadIdx.x < 2）再额外把 input[index + THREADS_PER_BLK] 写入 support[128 + threadIdx.x]（这是本 block 需要的 2 个跨块边界元素，即 halo）。总共执行 130 条 load 指令，而不是 v1 的 3×128 = 384 条。
步骤 2：__syncthreads() 保证所有 130 个数据都写进 shared memory 之后，任何线程才允许读取。
步骤 3：每个线程从 shared memory 读 3 个元素做累加（shared 是片上高速存储，比 DRAM 快得多）。
步骤 4：把结果写回 global memory。

【并行机制解说】

数据如何共享：support 是 per-block shared memory——所有线程用它交换数据，这正是”共享地址空间”编程风格（block 内）。
同步点：__syncthreads() 是唯一的块内屏障。它把”写入 shared”与”读取 shared”两个阶段隔开，防止线程 A 读到线程 B 尚未写入的旧数据。注意它只同步本 block 的线程，不同 block 之间没有任何同步语义（对应概念：__syncthreads、block 独立性）。
为什么省带宽：输入数据被复用了 3 次（相邻输出窗口重叠），v1 每次复用都走 DRAM，v2 只在开始时走一次 DRAM、之后全部命中片上 shared memory。这是”用 shared memory 显式管理局部性”的经典案例（对应概念：shared memory、memory bandwidth）。
运行在 V100 上：128 线程的 block 由 4 个 warp 执行（128/32）；由于块内线程必须”同时在场”（见陷阱 4），4 个 warp 必须同时驻留在同一 SM 上。
对应概念：shared memory、__syncthreads、SIMT/warp（4 warps/block）。

示例 4：用 shared memory 的树形归约（reduction）内核

// reduce.cu —— 编译：nvcc -o reduce reduce.cu
// 目标：把长度为 N 的数组 input 归约为一个总和（每个 block 产出一个部分和）
#include <cuda_runtime.h>

#define THREADS_PER_BLK 256
#define NUM_BLOCKS      64

__global__ void reduce(float* input, float* partial, int N)
{
    __shared__ float sdata[THREADS_PER_BLK];   // 每 block 一块共享累加区
    int tid = threadIdx.x;
    int i   = blockIdx.x * blockDim.x + tid;   // 该线程负责的全局元素下标

    // 1) 每线程把 global 中的一个元素装入 shared memory（越界元素当 0 处理）
    sdata[tid] = (i < N) ? input[i] : 0.f;
    __syncthreads();                           // 确保全部装入

    // 2) 树形归约：每轮参与线程数减半
    for (int s = THREADS_PER_BLK / 2; s > 0; s >>= 1) {
        if (tid < s)
            sdata[tid] += sdata[tid + s];      // 前半线程把后半的累加进来
        __syncthreads();                       // 防止读到自己还没被更新的数据
    }

    // 3) 每个 block 的 0 号线程把部分和写回 global
    if (tid == 0)
        partial[blockIdx.x] = sdata[0];
}

// host 端：把 NUM_BLOCKS 个部分和再相加（可再启动一个小 kernel，或拷回 CPU 求和）
int main()
{
    int N = NUM_BLOCKS * THREADS_PER_BLK;      // 64 × 256 = 16,384 个元素
    float *devInput, *devPartial;
    cudaMalloc(&devInput,    sizeof(float) * N);
    cudaMalloc(&devPartial,  sizeof(float) * NUM_BLOCKS);
    // （初始化 devInput 后……）
    reduce<<<NUM_BLOCKS, THREADS_PER_BLK>>>(devInput, devPartial, N);
    // （把 devPartial 拷回 host 求和，或再启动一个 reduce kernel 处理 64 个部分和）
    return 0;
}

【代码做了什么？】

步骤 1：256 个线程各装一个元素到 sdata，__syncthreads() 保证装载完成。
步骤 2：for (int s = 128; s > 0; s >>= 1) 执行 8 轮：第 1 轮线程 0-127 把 sdata[0..127] += sdata[128..255]，第 2 轮线程 0-63 把 sdata[0..63] += sdata[64..127]……每轮参与线程减半，8 轮后 sdata[0] 就是本 block 256 个元素的和。这是 O(log₂256) = 8 步的树形归约。
步骤 3：每个 block 由 0 号线程把部分和写入 partial[blockIdx.x]。
host 端最后再把 64 个部分和加起来（或再启动一次归约 kernel）。

【并行机制解说】

工作分配：每个 block 负责一块连续数据（256 个元素），block 之间完全独立（符合 CUDA 块可任意调度假设）。
数据共享与同步：sdata 是 shared memory；每轮加法后必须 __syncthreads()——否则线程 tid 可能读到相邻线程还没写完的旧值。归约是”写后读”依赖最密集的 kernel 之一，是理解 __syncthreads 的最佳练习。
性能观察：每轮有一半线程闲置（tid >= s 的线程直接跳到屏障）——这是朴素树归约的固有浪费，也是后续优化的方向（如 warp shuffle）。它还展示了”block 内 SPMD 协作”：256 个线程不是各自为战，而是作为一个整体协同计算（对应概念：block 内协作、shared memory、__syncthreads）。
注意：本 kernel 是教学版，实际还有”每线程加载多个元素（grid-stride loop）”等改进；它没有用原子操作，跨 block 的合并靠 host 完成。
对应概念：block 内 SPMD 协作、shared memory 归约、barrier。

示例 5（进阶/bonus）：persistent threads 编程风格

// persistent.cu —— 编译：nvcc -o persistent persistent.cu
// 思路：启动恰好"填满 GPU"的 block 数，让每个线程在 while 循环里反复取工作，
//       完全绕过 GPU 的线程块调度器，由应用自己管理工作分配。
#define THREADS_PER_BLK 128
#define BLOCKS_PER_CHIP 80 * (32 * 64 / 128)   // 针对 V100：80 SM × 2048 线程/SM ÷ 128

__device__ int workCounter = 0;                // global memory 中的工作计数器

__global__ void convolve(int N, float* input, float* output)
{
    __shared__ int startingIndex;              // 本 block 本轮起始下标（shared）
    __shared__ float support[THREADS_PER_BLK + 2];

    while (1) {
        if (threadIdx.x == 0)                  // 0 号线程原子地领取一段工作
            startingIndex = atomicInc(&workCounter, THREADS_PER_BLK);
        __syncthreads();                       // 广播给整个 block
        if (startingIndex >= N) break;         // 没有更多工作了

        int index = startingIndex + threadIdx.x;
        support[threadIdx.x] = input[index];
        if (threadIdx.x < 2)
            support[THREADS_PER_BLK + threadIdx.x] =
                input[index + THREADS_PER_BLK];
        __syncthreads();

        float result = 0.0f;
        for (int i = 0; i < 3; i++)
            result += support[threadIdx.x + i];
        output[index] = result;

        __syncthreads();                       // 防止下一轮覆盖 support 时还有人没读完
    }
}

// host：只启动 BLOCKS_PER_CHIP 个 block（不再 N/128 个）
// convolve<<<BLOCKS_PER_CHIP, THREADS_PER_BLK>>>(N, devInput, devOutput);

【代码做了什么？】

host 只启动恰好填满 GPU 的 block 数（V100 上 80 SM × 2048 线程/SM ÷ 128 = 1280 个 block）。
每个 block 的 0 号线程用 atomicInc 从全局计数器领取一段起始下标，__syncthreads() 广播后，整个 block 处理这 128 个输出；然后回到 while 循环领下一段，直到 startingIndex >= N。
工作分配从”硬件调度器”转移到”应用自身”——程序员的心理模型变成”所有 CUDA 线程同时跑在 GPU 上”。

【并行机制解说】

这是 CUDA 的反模式/特例：它要求程序员知道底层 GPU 的核数与每核容量（BLOCKS_PER_CHIP 写死了 V100 参数），并假设 GPU 确实会让所有 block 并发执行（”Ugg!”——幻灯片原话）。一旦换 GPU 或该假设不成立，程序性能甚至正确性都会出问题。
但它清楚地展示了两个概念：① block 内 __syncthreads() 充当”小循环内屏障”，把”领取工作→处理→再领取”串成流水；② 原子操作 atomicInc 是跨 block 共享 global 变量的唯一同步手段（对应概念：atomic、shared memory、block 调度）。
与”线程池”类比：就像 web server 启动时创建固定数量线程等待请求，而不是每来一个请求建一个线程——thread pool 的线程数是核数的函数，不是请求数的函数。

三、关键要点

CUDA 是”批量启动 + 两级层级”的数据并行编程：一次 kernel launch 创建成千上万个线程；问题分解成”block 集合”（网格），block 之间被假设无依赖、可按任意顺序调度（与 ISPC task 极其相似）；block 内部则是共享地址空间的 SPMD 编程（与 ISPC gang 相似）。没有 num_cores，程序天然可移植到不同规模 GPU。
warp/SIMT 是 CUDA 最重要的实现细节，但不是编程模型的一部分：32 个线程一个 warp，硬件运行时动态检查 warp 内指令是否一致，一致则用 SIMD ALU 一起执行；不一致则发散执行、性能受损。这与 ISPC 的”编译期生成 SIMD 指令”根本不同——CUDA 程序不会被编译成 SIMD 指令。
三种 device 地址空间 = 三种局部性：private（寄存器）/ shared（片上）/ global（DRAM）。用 shared memory 显式管理复用（如 1D 卷积把 384 次 global load 降到 130 次）是把 DRAM 带宽留给真正需要它的数据的关键。
block 内线程必须”同时在场”：因为 __syncthreads() 等机制允许块内线程互相依赖，系统不能先把 128 个线程跑完再跑后 128 个——block 开始时所有线程的寄存器上下文就必须全部分配好（这是对调度的硬约束，也是 block 大小受硬件上下文容量限制的原因）。
同步工具箱只有三样：__syncthreads()（块内屏障）、原子操作（global/shared）、kernel 返回的隐式全线程屏障（host/device 同步）。跨 block 的细粒度同步（如 spin-wait 等待别的 block）在 CUDA 里是危险甚至非法的，因为没有任何保证两个 block 会并发运行。

四、常见陷阱与注意事项

block 内线程数超过硬件上限：CUDA 限制每 block 最多 1024 个线程（V100 及多数现代 GPU）；shared memory 也是有限资源（V100 每 SM 128 KB shared+L1），block 占用过多资源会降低每 SM 可驻留 block 数，甚至无法启动（参见示例 3 的 fictitious core 只能放 2 个 block）。设计 block 大小时要同时考虑线程数、寄存器数、shared 字节数三个维度。
忘记 __syncthreads()（或放错位置）：shared memory 的”写后读”依赖必须用屏障隔开。漏掉屏障是 data race（读旧值）；把屏障放在条件分支里则可能导致死锁——因为屏障要求块内所有线程都到达，若只有部分线程执行了它，其余线程永远等不到。
shared memory 溢出 / 数组越界：support[THREADS_PER_BLK + 2] 这类”带 halo”的数组，越界写会悄悄破坏相邻 block 的数据（shared 是物理上连续的），且没有运行时错误提示——这是最难调的 bug 之一。
warp 发散（divergent execution）：if (threadIdx.x < 2) 这种分支会让 warp 内部分线程走 A 路径、部分走 B 路径，两条路径串行执行（掩码执行）。少量发散可接受，但大量发散（如按奇偶分叉）会让 SIMT 效率腰斩。设计时尽量让分支与 warp 边界对齐（如 if (warp_id == 0)）。
把 block 当 pthread 用：以为 block 之间可以像线程那样靠共享变量同步。记住 CUDA 只保证 block 可任意顺序调度——两个 block 用 spin-wait 互相等待（幻灯片中的 while(atomicAdd(&myFlag,0)==0){} 例子）在”每 SM 只驻留一个 block”的 GPU 上会永久死锁。跨 block 协作应通过 kernel 边界（启动多个 kernel）或原子操作完成。
host 端直接解引用 device 指针：devInput[i] 在 host 代码里是非法的（不同地址空间），必须 cudaMemcpy。反之亦然。
默认 launch 是同步的：host 调用 kernel 会阻塞直到 kernel 完成（CUDA 里其实是”默认串行流”语义），需要 host/device 并发时要用 streams——本讲没展开，但要知道 kernel 返回时刻有一个隐式全线程屏障。

五、思考题（带答案）

Q1：CUDA 是数据并行模型、共享地址空间模型，还是消息传递模型？请分别从”块之间”与”块内部”两个视角回答，并与 ISPC 的 gang/task 类比。

A1：答案在幻灯片总结里非常明确，是三者的”混合体”：

块之间（grid 层）：数据并行模型——问题被划分为独立 block，系统把它们调度到任意数量的核上（无 num_cores，类似 forall），block 间无依赖；这与 ISPC 的 task 几乎一模一样（task 也可按任意顺序调度）。
块内部（block 层）：共享地址空间模型的 SPMD 编程——线程并发运行、通过 shared memory 变量通信、用 __syncthreads() 同步；这与 ISPC 的 gang 类似，但 warp 不是编译期 SIMD（ISPC gang 编译成 SIMD 指令），而是硬件运行时动态检测的 SIMT。
host/device 之间：分布式地址空间——两个地址空间用 cudaMemcpy 搬数据，这又像消息传递模型。
所以：CUDA 一次编程同时用到了数据并行（块间）、共享内存（块内）、消息传递（host↔device）三种模型的元素——这正是幻灯片反复提醒”遇到任何并行编程系统先问它的语义是什么”的原因。

Q2：为什么 CUDA 必须为 block 内所有线程预先分配执行上下文（寄存器），而不能像 CPU 那样”先把 0-127 线程跑完再跑 128-255 线程”？

A2：因为 CUDA 允许（也鼓励）块内线程之间产生依赖——最简单的例子就是 __syncthreads()：如果线程 0-127 先跑完并越过了屏障，而 128-255 还没开始，那么”屏障”就名存实亡（前面的人不等后面的人）。更一般地，任何 shared memory 通信都假设双方同时存在。因此 CUDA 的语义是：block 开始执行时，块内所有线程都已存在且拥有寄存器状态；如果某线程可运行，它最终一定会被运行（不会死锁）。这给调度器套上了紧箍咒：一个 block 只能整体驻留在一个 SM 上，且占用的上下文数量受 SM 寄存器/线程槽容量限制——这也是为什么 block 不能无限大。

Q3：假设你要在 GPU 上统计一个数组的直方图（值域 0-9）。为什么 atomicAdd(&counts[A[i]], 1)（counts 在 global memory）是合法 CUDA 代码？它违反”块间无依赖”的假设吗？

A3：合法，且不违反块间无依赖假设。原子操作提供的是互斥（mutual exclusion）而不是顺序依赖——无论 block 以什么顺序、什么交错方式执行，每个 atomicAdd 都只做”读-加-写”的不可分割单元，最终 counts 的内容与执行顺序无关（都是把数组 A 中每个值出现次数统计一遍）。调度器仍然可以任意顺序调度 block。真正的麻烦是幻灯片中 block 0/block 1 用 while(atomicAdd(&myFlag,0)==0){} 互相 spin-wait 的代码：那是在要求特定的执行顺序（必须先跑 block 0），而 CUDA 不提供这种保证——在只放得下一个 block 的 GPU 上直接死锁。记住区分：原子操作 = 互斥（安全），顺序假设 = 依赖（危险）。

Q4（进阶）：为什么 warp 的指令流是”标量”指令，但执行效率却接近 SIMD？这和 ISPC 的 gang 有何区别？

A4：GPU 的每个硬件线程（warp 里的一个 lane）执行的是只有标量指令的指令流；硬件在取指时发现同一 warp 的 32 个线程正在执行同一条指令，就让 16 个 SIMD ALU 分两个时钟把它执行完（掩码掉不发散的 lane）。所以 SIMD 是运行时动态合并出来的，而不是编译期静态生成的。ISPC 则是编译器把 gang 内实例的操作静态编译成 SIMD 指令。区别的意义：CUDA 程序里写任意控制流（if、循环、函数调用）都不会”编译失败”，代价是发散时性能下降；而 ISPC 对 gang 内控制流的 SIMD 化是编译器的事。这也是为什么 CUDA 线程比 ISPC 实例”更自由”但需要程序员自觉保持 warp 一致性。

Share on

Bluesky Facebook LinkedIn Mastodon X (formerly Twitter)

Runguo Li