Vector Addition

agicy2026/6/6大约 5 分钟

Vector Addition

题目描述

编写一个 GPU 程序，对两个包含 32 位浮点数的向量执行逐元素加法。程序应接收两个长度相等的输入向量，并生成一个包含它们之和的输出向量。

实现要求

不允许使用外部库。
solve 函数签名必须保持不变。
最终结果必须存储在向量 $\mathbf{C}$ 中。

示例

示例 1

Input:  A = [1.0, 2.0, 3.0, 4.0]
        B = [5.0, 6.0, 7.0, 8.0]
Output: C = [6.0, 8.0, 10.0, 12.0]

示例 2

Input:  A = [1.5, 1.5, 1.5]
        B = [2.3, 2.3, 2.3]
Output: C = [3.8, 3.8, 3.8]

约束条件

输入向量 $\mathbf{A}$ 和 $\mathbf{B}$ 长度相同。
$1 \le N \le 100,000,000$ 。
性能测试将在 $N = 25,000,000$ 的规模下进行。

解题思路

向量加法是一个典型的数据并行（Data Parallelism）问题。由于每个元素的加法运算 $C[i] = A[i] + B[i]$ 都是独立的，我们可以利用 GPU 的大规模并行线程来同时计算多个元素。

模式分析

本题属于经典的 Map（映射）模式。在这种模式中，每个输入元素独立地通过相同的函数（此处为加法操作）映射到一个输出元素。由于计算过程没有数据依赖（如前缀和或归约），各线程之间完全不需要同步，这使得它成为最易于并行的算法结构。

从数据访问角度看，这是一个 Gather（聚集） 类型的操作：每个线程读取两个输入源 $A[i]$ 和 $B[i]$ ，并写入一个目标 $C[i]$ 。由于输入和输出索引一一对应且对齐，它也符合 Output-Centered（输出中心） 的设计思想，即以输出位置为锚点，分配线程去计算该位置的值。

效率分析

算术强度（Arithmetic Intensity）

向量加法的计算非常轻量，每个线程执行 1 次浮点加法操作，但需要进行 3 次 4 字节的全局内存访问（读取 $A[i]$ ，读取 $B[i]$ ，写入 $C[i]$ ）。

计算量：1 FLOP
访存量： $3 \times 4 = 12$ Bytes
算术强度： $1 / 12 \approx 0.083$ FLOPs/Byte

受限类型

由于算术强度极低，远低于任何现代 GPU 的峰值计算/带宽比（通常在 10~100+ FLOPs/Byte），因此这是一个典型的 Memory Bandwidth Bound（显存带宽受限） 内核。

程序的性能瓶颈完全取决于 GPU 的显存带宽，而非计算核心的算力。即使使用算力更强的 GPU，如果显存带宽没有提升，该程序的运行速度也不会有显著改善。优化方向应集中在确保内存合并访问（Coalesced Access）以最大化带宽利用率。

调度策略

采用 一维线性网格（1D Grid） 调度是最自然的选择。我们将 $N$ 个任务均匀分配给多个线程块（Thread Block）。

Grid 维度：根据向量总长度 $N$ 和每个 Block 的线程数动态计算。
Block 维度：通常选择 128、256 或 512 等 32 的倍数，以匹配 GPU 的 Warp 调度机制，最大化 Occupancy（占用率）。
线程索引：使用 $idx = blockIdx.x \times blockDim.x + threadIdx.x$ 计算全局唯一索引，直接对应数据下标。

这种简单的调度方式配合 GPU 的硬件调度器，能够自动处理成千上万个 Block 的执行顺序，无需用户手动管理任务队列。

核心策略

线程映射：将每个待计算的元素索引 $i$ 映射到一个唯一的 GPU 线程。
全局索引计算：利用 CUDA 内置变量计算当前线程对应的全局索引。

blockIdx.x：当前线程块的索引。
blockDim.x：每个线程块包含的线程数。
threadIdx.x：当前线程在块内的索引。
全局索引公式：idx = blockIdx.x * blockDim.x + threadIdx.x。

边界检查：确保计算出的索引 idx 小于向量长度 $N$ ，防止越界访问。

进阶优化：网格跨步循环 (Grid-Stride Loop)

虽然直接的一对一映射（每个线程处理一个元素）对于 $N$ 小于最大线程总数的情况是可行的，但在实际应用中， $N$ 往往远大于 GPU 能同时调度的线程数。此外，为了提高程序的健壮性和灵活性，推荐使用 网格跨步循环。

在内存带宽受限（Memory Bound）的场景下，通过限制 Block 的总数量并采用 Grid-Stride Loop，可以带来显著的资源优势：

减少上下文开销：GPU 需要为每个 Block 分配寄存器和共享内存资源。过多的 Block（尤其是当 $N$ 非常大时）会增加调度器的负担。限制 Grid 大小（例如设为 GPU SM 数量的 8~16 倍）可以减少这些管理开销。
避免资源争抢：在某些极端情况下，如果 Kernel 占用了过多的寄存器或共享内存，过大的 Grid 可能导致后续 Block 无法及时被调度，甚至在复杂的生产环境中与其他 Kernel 争抢资源。
提升指令级并行：每个线程串行处理多个元素，有助于掩盖指令流水线的延迟，增加单个线程的有效工作量。
原理：每个线程不仅处理索引 idx 的元素，还会在处理完当前元素后，跳过整个网格的大小（Grid Size），去处理下一个需要计算的元素 idx + gridDim.x * blockDim.x，直到覆盖所有数据。
优势：
解耦：代码逻辑与具体的 Grid 尺寸解耦，无论启动多少个 Block，程序都能正确运行。
重用：增加了每个线程的工作量，有助于掩盖指令延迟和初始化开销。
调试：即使配置为 1 个 Block 1 个 Thread，程序也能串行正确执行，便于调试。

代码实现

CUDA

#include <cuda_runtime.h>

__global__ void vector_add(const float* A, const float* B, float* C, int N) {
    int i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < N) C[i] = A[i] + B[i];
}

// Grid-Stride Loop
__global__ void vector_add_grid_stride(const float* A, const float* B, float* C, int N) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    int stride = gridDim.x * blockDim.x;
    for (int i = idx; i < N; i += stride) C[i] = A[i] + B[i];
}

// float4 vectorized (128-bit load/store)
__global__ void vector_add_float4(const float* A, const float* B, float* C, int N) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    int stride = gridDim.x * blockDim.x;
    int N4 = N / 4;
    const float4* A4 = (const float4*)A;
    const float4* B4 = (const float4*)B;
    float4* C4 = (float4*)C;
    for (int i = idx; i < N4; i += stride) {
        float4 a = A4[i], b = B4[i];
        float4 c = {a.x+b.x, a.y+b.y, a.z+b.z, a.w+b.w};
        C4[i] = c;
    }
    for (int i = N4*4 + idx; i < N; i += stride) C[i] = A[i] + B[i];
}

extern "C" void solve(const float* A, const float* B, float* C, int N) {
    int threadsPerBlock = 256;
    int blocksPerGrid = (N + threadsPerBlock - 1) / threadsPerBlock;
    vector_add<<<blocksPerGrid, threadsPerBlock>>>(A, B, C, N);
    cudaDeviceSynchronize();
}

Triton

import triton, triton.language as tl, torch

@triton.jit
def vector_add_kernel(A_ptr, B_ptr, C_ptr, N: tl.constexpr, BLOCK_SIZE: tl.constexpr):
    idx = tl.program_id(0) * BLOCK_SIZE + tl.arange(0, BLOCK_SIZE)
    mask = idx < N
    a = tl.load(A_ptr + idx, mask=mask)
    b = tl.load(B_ptr + idx, mask=mask)
    tl.store(C_ptr + idx, a + b, mask=mask)

def solve(A, B):
    N = A.numel()
    C = torch.empty_like(A)
    grid = lambda meta: (triton.cdiv(N, meta["BLOCK_SIZE"]),)
    vector_add_kernel[grid](A, B, C, N, BLOCK_SIZE=1024)
    return C