深入理解CUDA：统一计算设备架构

CUDA（Compute Unified Device Architecture，统一计算设备架构）是NVIDIA开发的一个并行计算平台和编程模型，于2006年推出。它允许开发者使用C/C++、Fortran、Python等编程语言来利用NVIDIA GPU的并行计算能力。

1. CUDA架构概述

1.1 设计理念

CUDA的设计核心理念是提供一种通用的并行计算架构，让开发者能够：

高效利用GPU的大规模并行计算能力：现代GPU拥有数千个计算核心，能够同时执行数万个线程
保持编程的易用性：通过扩展C/C++语言，使开发者能够以熟悉的方式编写GPU程序
实现CPU与GPU的协同工作：CUDA提供了heterogeneous programming（异构编程）模型

1.2 硬件架构

现代NVIDIA GPU由以下关键组件组成：

GPU
├── Streaming Multiprocessors (SMs)
│   ├── CUDA Cores
│   ├── Special Function Units (SFUs)
│   ├── Load/Store Units
│   └── Shared Memory
├── Memory Controllers
├── Global Memory (DRAM)
└── L2 Cache

1.3 计算能力（Compute Capability）

每个GPU都有一个计算能力版本号，表示其支持的硬件特性：

Compute Capability 3.x: Kepler架构，引入动态并行
Compute Capability 5.x: Maxwell架构，改进内存管理
Compute Capability 6.x: Pascal架构，引入统一内存和NVLink
Compute Capability 7.x: Volta和Turing架构，引入Tensor Cores
Compute Capability 8.x: Ampere架构，增强Tensor Cores性能
Compute Capability 9.x: Hopper架构，专为HPC和AI优化

2. CUDA编程模型

2.1 核心概念

Kernel函数

Kernel是在GPU上执行的函数，由大量线程并行执行：

__global__ void kernelName(float* input, float* output, int n) {
    int tid = blockIdx.x * blockDim.x + threadIdx.x;
    if (tid < n) {
        output[tid] = input[tid] * 2.0f;
    }
}

线程层次结构（Thread Hierarchy）

CUDA采用三级线程层次结构：

Grid (网格)
├── Block 0 (线程块)
│   ├── Thread 0-1023
│   └── ...
├── Block 1
│   └── ...
└── Block N

Thread（线程）: 最基本的执行单元
Block（线程块）: 一组可以协作的线程（最多1024个）
Grid（网格）: 由多个线程块组成的完整kernel执行

执行配置

使用<<<blocks_per_grid, threads_per_block>>>语法启动kernel：

// 启动256个线程块，每个线程块512个线程
kernelName<<<256, 512>>>(d_input, d_output, n);

2.2 内置变量

CUDA提供以下内置变量用于定位线程：

// 线程索引
threadIdx.x, threadIdx.y, threadIdx.z
// 块索引
blockIdx.x, blockIdx.y, blockIdx.z
// 块维度
blockDim.x, blockDim.y, blockDim.z
// 网格维度
gridDim.x, gridDim.y, gridDim.z

// 计算全局线程ID
int globalId = blockIdx.x * blockDim.x + threadIdx.x;

3. 内存层次结构

3.1 内存类型与特性

CUDA提供多种内存类型，具有不同的访问速度和作用域：

寄存器（Registers）

位置: 芯片内，每个线程私有
访问速度: 最快
容量: 有限（每个线程约64个寄存器）
生命周期: 线程生命周期内有效

__global__ void kernel() {
    int temp = 10; // 寄存器变量
    // ...
}

局部内存（Local Memory）

位置: 全局内存中，每个线程私有
访问速度: 慢（需要通过内存控制器访问）
用途: 当寄存器不足时，编译器自动将变量溢出到局部内存

共享内存（Shared Memory）

位置: 芯片内，线程块内共享
访问速度: 快（接近寄存器速度）
容量: 配置有限（每个SM约48KB）
生命周期: 线程块生命周期内有效

__global__ void kernel() {
    __shared__ float sdata[256]; // 共享内存声明
    // 线程块内协作处理
    sdata[threadIdx.x] = global_data[blockIdx.x * blockDim.x + threadIdx.x];
    __syncthreads(); // 同步线程块内所有线程
}

全局内存（Global Memory）

位置: 板载DRAM，所有线程可访问
访问速度: 较慢（高延迟，但高带宽）
容量: 大（GB级别）
持久性: 整个应用程序生命周期内有效

// 在主机端分配全局内存
float *d_data;
cudaMalloc(&d_data, size * sizeof(float));

// 访问全局内存
__global__ void kernel(float* d_data) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    float value = d_data[idx]; // 从全局内存读取
    // ...
}

常量内存（Constant Memory）

位置: 片内缓存，只读
访问速度: 快（缓存命中时），广播优化
容量: 64KB
用途: 存放不变的数据，内核参数等

__constant__ float const_data[1024]; // 常量内存声明

// 在主机端设置常量内存
cudaMemcpyToSymbol(const_data, h_data, size);

纹理内存（Texture Memory）

位置: 带缓存的2D/3D内存空间
访问速度: 快（空间局部性优化）
用途: 图像处理、空间数据访问模式

3.2 内存优化策略

合并访问（Coalesced Access）

合并访问是CUDA性能优化的关键技术，指warp内线程连续访问全局内存：

// 良好的合并访问
__global__ void goodCoalesced(float* input, float* output, int n) {
    int tid = blockIdx.x * blockDim.x + threadIdx.x;
    output[tid] = input[tid] * 2.0f;
}

// 不好的非合并访问
__global__ void badCoalesced(float* input, float* output, int n) {
    int tid = blockIdx.x * blockDim.x + threadIdx.x;
    output[tid] = input[tid * 16]; // 跳跃访问
}

共享内存利用

使用共享内存减少全局内存访问：

__global__ void matrixMultiply(float* A, float* B, float* C, int N) {
    __shared__ float As[BLOCK_SIZE][BLOCK_SIZE];
    __shared__ float Bs[BLOCK_SIZE][BLOCK_SIZE];

    int bx = blockIdx.x, by = blockIdx.y;
    int tx = threadIdx.x, ty = threadIdx.y;

    // 加载块到共享内存
    int row = by * BLOCK_SIZE + ty;
    int col = bx * BLOCK_SIZE + tx;

    for (int k = 0; k < N / BLOCK_SIZE; k++) {
        As[ty][tx] = A[row * N + k * BLOCK_SIZE + tx];
        Bs[ty][tx] = B[(k * BLOCK_SIZE + ty) * N + col];
        __syncthreads();

        // 计算乘积
        float sum = 0;
        for (int i = 0; i < BLOCK_SIZE; i++) {
            sum += As[ty][i] * Bs[i][tx];
        }
        C[row * N + col] = sum;
        __syncthreads();
    }
}

4. CUDA编程范式

4.1 异构编程模型

CUDA程序的基本结构：

#include <cuda_runtime.h>
#include <stdio.h>

// 核函数（在GPU上执行）
__global__ void add(float* a, float* b, float* c, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < n) {
        c[idx] = a[idx] + b[idx];
    }
}

int main() {
    int n = 1024;
    size_t size = n * sizeof(float);

    // 主机端数据
    float *h_a = (float*)malloc(size);
    float *h_b = (float*)malloc(size);
    float *h_c = (float*)malloc(size);

    // 设备端数据
    float *d_a, *d_b, *d_c;
    cudaMalloc(&d_a, size);
    cudaMalloc(&d_b, size);
    cudaMalloc(&d_c, size);

    // 数据传输：主机到设备
    cudaMemcpy(d_a, h_a, size, cudaMemcpyHostToDevice);
    cudaMemcpy(d_b, h_b, size, cudaMemcpyHostToDevice);

    // 启动核函数
    int blockSize = 256;
    int gridSize = (n + blockSize - 1) / blockSize;
    add<<<gridSize, blockSize>>>(d_a, d_b, d_c, n);

    // 数据传输：设备到主机
    cudaMemcpy(h_c, d_c, size, cudaMemcpyDeviceToHost);

    // 清理资源
    cudaFree(d_a);
    cudaFree(d_b);
    cudaFree(d_c);
    free(h_a);
    free(h_b);
    free(h_c);

    return 0;
}

4.2 错误处理

CUDA错误处理最佳实践：

#define CUDA_CHECK(call) \
    do { \
        cudaError_t error = call; \
        if (error != cudaSuccess) { \
            fprintf(stderr, "CUDA error at %s:%d - %s\n", \
                    __FILE__, __LINE__, cudaGetErrorString(error)); \
            exit(EXIT_FAILURE); \
        } \
    } while (0)

// 使用示例
CUDA_CHECK(cudaMalloc(&d_ptr, size));
CUDA_CHECK(cudaMemcpy(d_ptr, h_ptr, size, cudaMemcpyHostToDevice));
kernel<<<grid, block>>>(d_ptr);
CUDA_CHECK(cudaGetLastError()); // 检查kernel启动错误
CUDA_CHECK(cudaDeviceSynchronize()); // 等待kernel完成

4.3 流（Streams）与并发

CUDA流支持异步操作和并发执行：

// 创建多个流
cudaStream_t stream1, stream2;
cudaStreamCreate(&stream1);
cudaStreamCreate(&stream2);

// 异步内存拷贝和kernel执行
cudaMemcpyAsync(d_a1, h_a1, size, cudaMemcpyHostToDevice, stream1);
cudaMemcpyAsync(d_a2, h_a2, size, cudaMemcpyHostToDevice, stream2);

kernel<<<grid1, block1, 0, stream1>>>(d_a1, d_b1, d_c1);
kernel<<<grid2, block2, 0, stream2>>>(d_a2, d_b2, d_c2);

cudaMemcpyAsync(h_c1, d_c1, size, cudaMemcpyDeviceToHost, stream1);
cudaMemcpyAsync(h_c2, d_c2, size, cudaMemcpyDeviceToHost, stream2);

// 同步所有流
cudaStreamSynchronize(stream1);
cudaStreamSynchronize(stream2);

// 清理资源
cudaStreamDestroy(stream1);
cudaStreamDestroy(stream2);

5. 性能优化技术

5.1 占用率优化（Occupancy）

占用率是衡量SM上活动warp数量的指标：

// 查询设备属性
cudaDeviceProp prop;
cudaGetDeviceProperties(&prop, 0);

// 计算最佳block size
int minGridSize, blockSize;
cudaOccupancyMaxPotentialBlockSize(&minGridSize, &blockSize,
                                  myKernel, 0, 0);

// 启动优化后的kernel
myKernel<<<minGridSize, blockSize>>>(/* parameters */);

5.2 指令级优化

循环展开

// 展开前
__global__ void unroll1(float* data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    for (int i = 0; i < 4; i++) {
        data[idx * 4 + i] *= 2.0f;
    }
}

// 展开后
__global__ void unroll2(float* data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    int base = idx * 4;
    data[base] *= 2.0f;
    data[base + 1] *= 2.0f;
    data[base + 2] *= 2.0f;
    data[base + 3] *= 2.0f;
}

使用限制器和predication

// 使用min避免分支
__global__ void optimized(float* data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    int stride = blockDim.x * gridDim.x;

    for (int i = idx; i < n; i += stride) {
        data[i] = sqrt(data[i]);
    }
}

5.3 内存访问优化

矢量化访问

// 使用float4进行矢量化访问
__global__ void vectorizedLoad(float4* input, float4* output, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < n) {
        float4 data = input[idx];
        data.x *= 2.0f;
        data.y *= 2.0f;
        data.z *= 2.0f;
        data.w *= 2.0f;
        output[idx] = data;
    }
}

只读数据缓存

// 使用__ldg读取只读数据
__global__ void readOnlyCache(const float* __restrict__ input,
                             float* output, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < n) {
        float val = __ldg(&input[idx]); // 使用只读缓存
        output[idx] = val * 2.0f;
    }
}

6. 高级特性

6.1 动态并行（Dynamic Parallelism）

允许GPU上运行的kernel启动其他kernel：

// 子kernel
__global__ void childKernel(float* data, int offset, int size) {
    int idx = threadIdx.x;
    if (idx < size) {
        data[offset + idx] *= 2.0f;
    }
}

// 父kernel
__global__ void parentKernel(float* data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;

    // 动态启动子kernel
    childKernel<<<1, 256>>>(data, idx * 256, 256);
    cudaDeviceSynchronize(); // 同步子kernel
}

6.2 统一内存（Unified Memory）

简化内存管理，自动处理数据迁移：

// 统一内存分配
float *unified_data;
cudaMallocManaged(&unified_data, size * sizeof(float));

// 直接访问，无需显式拷贝
kernel<<<grid, block>>>(unified_data, size);

// CPU也可以直接访问
for (int i = 0; i < size; i++) {
    unified_data[i] = i * 2.0f;
}

cudaFree(unified_data);

6.3 图（Graphs）

优化重复任务的启动开销：

cudaGraph_t graph;
cudaGraphExec_t graphExec;

// 创建图
cudaStreamBeginCapture(stream, cudaStreamCaptureModeGlobal);
kernel<<<grid1, block1, 0, stream>>>(/* params */);
cudaMemcpyAsync(d_ptr, h_ptr, size, cudaMemcpyHostToDevice, stream);
kernel<<<grid2, block2, 0, stream>>>(/* params */);
cudaStreamEndCapture(stream, &graph);

// 实例化图
cudaGraphInstantiate(&graphExec, graph, NULL, NULL, 0);

// 重复执行图
for (int i = 0; i < iterations; i++) {
    cudaGraphLaunch(graphExec, stream);
    cudaStreamSynchronize(stream);
}

7. 性能分析工具

7.1 Nsight Compute

细粒度kernel性能分析：

# 分析kernel性能
nv-nsight-cu ./my_application

# 指定特定kernel
nv-nsight-cu --kernel=my_kernel ./my_application

7.2 Nsight Systems

系统级性能分析：

# 生成系统级跟踪
nsys profile ./my_application

# 查看结果
nsys-ui profile.nsys-rep

7.3 CUDA Profiler

使用CUDA profiler API：

cudaProfilerStart();
kernel<<<grid, block>>>(/* parameters */);
cudaProfilerStop();

8. 最佳实践

8.1 通用优化原则

最大化并行度：确保有足够的活跃线程
优化内存访问：合并访问、使用共享内存
减少分支发散：避免warp内线程执行不同路径
平衡计算与内存访问：避免内存瓶颈

8.2 代码组织

// 检查编译时常量
#if !defined(__CUDA_ARCH__) || __CUDA_ARCH__ >= 600
    // 使用新特性
#else
    // 兼容旧架构
#endif

// 模板化不同数据类型
template<typename T>
__global__ void genericKernel(T* data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < n) {
        data[idx] = process(data[idx]);
    }
}

8.3 调试技巧

// 编译时启用调试
nvcc -g -G my_code.cu

// 使用printf调试（谨慎使用，影响性能）
__global__ void debugKernel(float* data) {
    printf("Thread %d: data = %f\n", threadIdx.x, data[threadIdx.x]);
}

// 使用assert进行运行时检查
__global__ void kernelWithAssert(int* data) {
    int idx = threadIdx.x;
    assert(idx < 1024);
    data[idx] = idx;
}

9. 应用领域

CUDA在以下领域得到广泛应用：

9.1 科学计算

流体动力学：CFD模拟
分子动力学：蛋白质折叠模拟
天体物理：N体问题求解

9.2 人工智能与深度学习

神经网络训练：CNN、RNN、Transformer
推理加速：实时AI应用
强化学习：大规模策略搜索

9.3 图像与视频处理

计算机视觉：目标检测、图像分割
视频编解码：实时处理
图像渲染：光线追踪

9.4 金融计算

风险管理：蒙特卡洛模拟
期权定价：Black-Scholes模型
高频交易：算法优化

10. 总结

CUDA提供了一个强大而灵活的并行计算平台，通过以下关键特性实现了GPU计算能力的高效利用：

层次化的编程模型：Grid-Block-Thread结构提供了灵活的并行组织方式
多级内存层次：不同速度和容量的内存类型优化了数据访问
丰富的API和工具：从基础API到高级库，从编译器到性能分析工具
持续的架构演进：从Tesla到Hopper，不断引入新特性提升性能

成功的CUDA编程需要深入理解GPU架构特性，合理设计并行算法，并进行系统性的性能优化。随着NVIDIA GPU架构的不断演进，CUDA将继续在HPC和AI领域发挥重要作用。