引言
作为一名使用eBPF进行CPU性能分析的工程师，在转向学习GPU性能优化分析时，一直在思考GPU上是否有技术也可以实现用户自定义探针式性能分析？学习NVIDIA Nsight Compute过程中我发现，尽管它提供了丰富的硬件计数器和细粒度的性能指标，但是其数据采集范围和触发条件在 profiling 前由硬件或Metrics固定，用户可控程度有限。Nsight Compute虽然会给我一个非常全局的Kernel性能报告，但是学习过程中我会疑问，如果我想知道特定条件下更细粒度的性能数据，比如某个Kernel执行时第 128 个warpid到底花了多少时间，好像它没法告诉我？
最近一个基于PTX的动态插桩框架（Neutrino）引起了我的注意，它允许我们在GPU Kernel编译后的PTX代码中，动态插入探针指令，在运行时收集自定义指标，比如特定变量的值、内存访问模式、耗时统计等，这是一种可编程的、条件驱动的观测工具，可以满足我的细粒度、个性化的数据捕获需求，在PTX层实现了“GPU版eBPF”。
虽然发现它还无法像eBPF那样完全动态attach到正在运行的Kernel，但已经实现了“编译后插桩 + 运行时观测”的核心能力，这让我非常有兴趣研究研究它，本文以总结自己学习成果为契机，带着各位看官一起来认识GPU Kernel的性能分析世界，从 CUDA 基础讲起，通过一个粗糙的矩阵乘法例子，逐步展示如何结合Nsight Compute、PTX插桩与Neutrino的分析框架，探下GPU Kernel性能分析的深度。
本人刚开始入门学习GPU性能调优相关知识，所以若文章中有不够深刻或者理解偏差的部分，欢迎大家留言～
第一部分：CUDA架构基础，是理解性能瓶颈的起点
万丈高楼平地起，展开主线内容前，想先聊聊CUDA基础，因为它们是识别与解决性能瓶颈的理论依据。
GPU的高性能源于其大规模并行处理能力，要进行有效的性能分析得先掌握其硬件架构和编程模型。
1. GPU核心结构：SM、warp、thread/block/grid 的调度关系
左侧是逻辑实现，右侧是GPU硬件架构，由下往上，从左到右地分层说明：
1. Grid代表Kernel的单次完整调用，当一次Kernel被调用时，其所有线程（整个Grid）都由单一GPU设备负责执行
2. Grid里包含若干Block（线程块），这些Block会被分配给设备上可用的SM（流式多处理器），一个Block内的所有Thread（线程）被保证可以调度到同一个SM上，不可拆分到其他SM上执行。SM会将Block中的线程重新分组成大小固定，通常为32个线程的束(Warp)，束是GPU硬件调度的基本单位，SM上的Warp调度器以Warp为单位进行指令的发送和管理，一个 Warp 内的 32个线程在理想情况下，同一周期会执行相同的指令，这种执行模式被称为SIMT (Single Instruction, Multiple Threads)
3. Thread是最基本的执行单位，由SM里的一个个计算逻辑单元（ALU）执行。每一个Thread都会把Kernel从头到尾完整地执行一遍，每个Thread各自拥有完全独立的Registers（寄存器），并且这些寄存器里保存着各自不同的值。
2. 内存层级：Global、Shared、Registers、Cache的区别
如上图展示，CUDA架构定义了一个分层的内存模型，不同层级的内存具有不同的访问速度、容量和作用域，这是性能优化的关键考量因素：
3. 通过具体示例理解概念与性能的直接关联
Kernel性能调优围绕四个基本策略展开：
- 最大化并行执行以实现最大利用率
- 多kernel并发执行+CUDA 流，最大化GPU 各个SM利用率
- 线程级并行，最大化ALU利用率，通过Warp调度实现延迟隐藏
- 指令级并行，减少指令间依赖性，隐藏指令延迟
- 优化内存使用率以实现最大内存吞吐量
- 最小化主机内存与设备内存间的数据传输
- 最大限度地利用片上内存来最小化全局内存使用
- 通过合并内存访问，提高全局内存吞吐
- 共享内存访问减少bank conflict
- 优化指令使用以实现最大指令吞吐量
- 减少吞吐量低的算术指令的使用
- 避免线程束发散
- 最小化内存抖动，避免频繁的分配和释放内存
因为以上优化策略完全展开的话篇幅过长，所以只是点到为止，本小节只挑选部分内容具体介绍，对以上优化策略感兴趣的看官可以阅读 CUDA C++编程指南。
典型示例
这是一个使用CUDA C++实现的矩阵乘法例子，之所以挑选矩阵乘法作为示例，是因为这在深度学习中最常见，它既有计算，又有内存访问，非常适合用来观察性能瓶颈。
#include<iostream> #include<vector> #include<random> #include<stdexcept> #include<cmath> #include<cuda_runtime.h> #define CUDA_CHECK(err) { \ cudaError_t err_ = (err); \ if (err_ != cudaSuccess) { \ std::cerr << "CUDA error at " << __FILE__ << ":" << __LINE__; \ std::cerr << ": " << cudaGetErrorString(err_) << std::endl; \ exit(EXIT_FAILURE); \ } \ } voidgemm_cpu(constfloat* A, constfloat* B, float* C, int M, int N, int K){ for (int r = 0; r < M; ++r) { for (int c = 0; c < N; ++c) { float sum = 0.0f; for (int i = 0; i < K; ++i) { sum += A[r * K + i] * B[i * N + c]; } C[r * N + c] = sum; } } } voidverify_result(constfloat* C_gpu, constfloat* C_cpu, int M, int N){ double max_error = 0.0; for (int i = 0; i < M * N; ++i) { max_error = fmax(max_error, fabs(C_gpu[i] - C_cpu[i])); } std::cout << "Max error between GPU and CPU results: " << max_error << std::endl; if (max_error > 1e-4) { std::cout << "Verification FAILED!" << std::endl; } else { std::cout << "Verification PASSED!" << std::endl; } } __global__ voidgemm_gpu(constfloat* A, constfloat* B, float* C, int M, int N, int K){ int row = blockIdx.y * blockDim.y + threadIdx.y; int col = blockIdx.x * blockDim.x + threadIdx.x; if (row < M && col < N) { float sum = 0.0f; for (int i = 0; i < K; ++i) { sum += A[row * K + i] * B[i * N + col]; } C[row * N + col] = sum; } } intmain(int argc, char** argv){ int M = 1024; int N = 1024; int K = 1024; std::cout << "Matrix dimensions: C(" << M << "x" << N << ") = A(" << M << "x" << K << ") * B(" << K << "x" << N << ")" << std::endl; std::vector<float> h_A(M * K); std::vector<float> h_B(K * N); std::vector<float> h_C_gpu(M * N, 0); std::vector<float> h_C_cpu(M * N, 0); std::random_device rd; std::mt19937 gen(rd()); std::uniform_real_distribution<> dis(0.0, 1.0); for (int i = 0; i < M * K; ++i) h_A[i] = dis(gen); for (int i = 0; i < K * N; ++i) h_B[i] = dis(gen); float *d_A, *d_B, *d_C; CUDA_CHECK(cudaMalloc(&d_A, M * K * sizeof(float))); CUDA_CHECK(cudaMalloc(&d_B, K * N * sizeof(float))); CUDA_CHECK(cudaMalloc(&d_C, M * N * sizeof(float))); CUDA_CHECK(cudaMemcpy(d_A, h_A.data(), M * K * sizeof(float), cudaMemcpyHostToDevice)); CUDA_CHECK(cudaMemcpy(d_B, h_B.data(), K * N * sizeof(float), cudaMemcpyHostToDevice)); int recommend_blocksize; int recommend_gridsize; cudaOccupancyMaxPotentialBlockSize( &recommend_gridsize, &recommend_blocksize, (void*)gemm_gpu, 0, 0); int device; cudaDeviceProp prop; cudaGetDevice(&device); cudaGetDeviceProperties(&prop, device); int dimX = prop.warpSize; int dimY = 1; if (recommend_blocksize % dimX == 0){ dimY = recommend_blocksize/dimX; } else{ dimX = static_cast<int>(sqrt(recommend_blocksize)); dimY = static_cast<int>(sqrt(recommend_blocksize)); } dim3 block_dim(dimX, dimY); dim3 grid_dim((N + block_dim.x - 1) / block_dim.x, (M + block_dim.y - 1) / block_dim.y); std::cout << "\n--- Running GEMM ---" << std::endl; gemm_gpu<<<grid_dim, block_dim>>>(d_A, d_B, d_C, M, N, K); CUDA_CHECK(cudaDeviceSynchronize()); CUDA_CHECK(cudaMemcpy(h_C_gpu.data(), d_C, M * N * sizeof(float), cudaMemcpyDeviceToHost)); gemm_cpu(h_A.data(), h_B.data(), h_C_cpu.data(), M, N, K); verify_result(h_C_gpu.data(), h_C_cpu.data(), M, N); CUDA_CHECK(cudaFree(d_A)); CUDA_CHECK(cudaFree(d_B)); CUDA_CHECK(cudaFree(d_C)); return0; }
代码结构简析
代码实现的是一个GPU加速的矩阵乘法计算，并通过与CPU计算结果对比，验证正确性。计算公式为：C = A × B，其中：
- A 是 M×K 矩阵
- B 是 K×N 矩阵
- C 是 M×N 结果矩阵
开始 │ ├─ line55-57 --定义矩阵大小 (M, N, K) │ ├─ line61-70 --主机内存分配与初始化 │ ├─ line72-75 --GPU 内存分配 (d_A, d_B, d_C) │ ├─ line76-77 --数据从主机传送到设备 (H→D) │ ├─ line79-101 --评估最佳block_size，配置block_dim、grid_dim │ ├─ line104 --启动Kernel gemm_gpu <<<grid, block>>> │ ├─ line42-52 --GPU上Kernel执行 │ ├─ line105 --GPU同步等待完成 │ ├─ line106 --结果从设备传回主机 (D→H) → h_C_gpu │ ├─ line107、17-27 --CPU计算参考 gemm_cpu() │ ├─ line108、29-40 --对比CPU、GPU计算结果 verify_result() │ ├─ line110-112 --释放GPU内存 │ 程序结束
CUDA编程模型核心
这里是结合示例代码简要概括CUDA的编程模型要点，理解它对后续章节内容也是必要的：
- 算子是通过
kernel_name<<>>(argument list)方式启动的，其中grid表示有多少个block，block表示一个block里有多少个thread，从而可知这个kernel使用的线程总数。GPU并行化思路是让每一个线程负责矩阵C中的一个元素的计算，比如线程mn处理C(m,n)，因此线程总数与M*N的矩阵元素个数进行适配，每个线程都会实际的执行一遍gemm_gpu里的代码。 - 以我的实验环境为例，采用CUDA计算能力为8.6的GPU卡，每SM最大线程数1536，若block_dim（16,16）256个线程，一个SM上可以运行6个Block，每个Block占用8个Warp。若block_dim（32,32）1024个线程，由于单SM最大线程数的限制，一个SM上就只能运行1个Block。示例通过
cudaOccupancyMaxPotentialBlockSize评估最佳block_dim大小，规划对于C矩阵的线程布局。 - 设备内存是线性的，每个线程具备唯一标识，假设一个以<<<4096,256>>>启动的kernel，线程布局如图所示，通过
blockIdx.i * blockDim.i + threadIdx.i定位到具体的线程ID，每一个线程将会去计算一个C矩阵中的一个元素的值。
- gemm_gpu函数体for循环内部则是具体矩阵乘法的点积过程。如图所示，2*3的矩阵A与3*2的矩阵B相乘，对于得到的2*2的矩阵C来说，C的第m行第n列元素是A的第m行每个元素与矩阵B到第n列对应的每个元素的乘积之和。抽象为数学公式的话，对于A=M*K，B=K*N，C=M*N而言，C[m][n] = Σ (A[m][i] * B[i][n])，i为从0->K-1，第一次是A的第m行第一个元素乘B的第一行第n列元素，第一次是A的第m行第二个元素乘B的第二行第n列元素，以此类推
以上是矩阵点积的数学表达，二维矩阵在扁平化的一维内存中的寻址方式是按照行主序排列：先存放第0行的所有元素，然后是第1行的所有元素，接着是第2行，以此类推。所以对于上图2*3的A矩阵，在内存中被“铺平”成下图这样。
内存地址 -> [ a00, a01, a02, a10, a11, a12 ] ^---------------^ ^----------- ^ 第0行 (3个元素) 第1行 (3个元素)
对于具体的一个线程在计算给定的C(row,col)时，过程如下图：
性能调优
以上四种性能优化策略方向，具体哪种策略能够为程序带来最佳性能提升，取决于程序性能受限原因：如果对性能瓶颈在内存的kernel的指令一个劲优化，并不会带来显著的性能提升，因为性能瓶颈不在指令层面。
把示例代码编译为可执行程序，然后通过Nsight Compute分析出性能可改进点最为高效。
经分析该kernel最大的性能瓶颈是内存访问效率低（非合并访问），这浪费了大量的内存带宽，有16.7%的理论优化空间，其次内存指令队列有拥堵现象（全局内存访问频繁），导致了轻微的停顿：
减少非合并访问与全局内存访问频次
非合并访问
Block里的Thread在Warp上的分布是有序的，通常WarpSize是32，每个Warp包含32个连续线程ID的线程。合并访问指的是Warp级的合并内存访问，在处理同一条指令下，一个Warp的32个线程同时访问32个连续地址，就能将访问请求合并为少数次的内存事务统一提交给硬件处理，由此提升内存吞吐。
下图这是举例BlockDim（16，16）里的每个线程在Warp的分布，用此来反观示例代码的内存寻址为什么是非合并访问：
threadIdx.x-> 0 1 2 3 ... 15 +---+---+---+---+---+ 0 | 0 | 1 | 2 | 3 |...| 15 ← 第0行 → 线性ID: 0~15 → warp 0 ty +---+---+---+---+---+ 1 |16 |17 |18 |19 |...| 31 ← 第1行 → 线性ID: 16~31 → warp 0 +---+---+---+---+---+ 2 |32 |33 |... | ← 第2行 → 线性ID: 32~47 → warp 1 +---+---+---+---+---+ 3 |48 |... | ← 第3行 → 线性ID: 48~63 → warp 1 +---+---+---+---+---+ 4 |64 |... | ← 第4行 → 线性ID: 64~79 → warp 2 +---+---+---+---+---+ ...
1.对于A[row * K + i]，它的地址由row决定，row = blockIdx.y * blockDim.y + threadIdx.y，对warp0内线程来说，前16个线程 (threadIdx.y =0)，它们的row值是完全相同的，后16个线程 (threadIdx.y =1)，它们的row值也是完全相同的。其中线程ID15、16，这里跨度了K*4字节地址空间，这造成4M的跨度，硬件无法将这些访问合并；并且对前16个线程，每个线程地址实际访问地址一样，但都需要独立从全局内存加载数据，造成内存浪费，同理后16个线程也存在内存浪费。
2.对于B[i * N + col]，它的地址由col决定，col = blockIdx.x * blockDim.x + threadIdx.x，对warp0内前16个线程，threadIdx.x是连续的因此这16个线程是访问的连续的16*4字节地址空间，同理后16个线程也将访问连续的16*4字节地址空间，满足内存合并要求。
对于示例中BlockDim（32，24）而言同理，A矩阵的寻址是非合并访问，是需要进行优化的地方。
全局内存访问时延
这是一份公开的基于A100的GPU microbenchmark基准测试，与我实验环境GPU型号相近，可见内存访问延迟呈现出巨大的层级差异。
全局内存访问约290个指令时钟周期-CPI（cycles per instruction），对于示例中一次add.s64指令花费4个CPI来说，一次对全局内存的访问所花时间是一次指令执行的72.5倍。
当一个Kernel中的所有线程都需要从全局内存加载数据时，Warp调度器将没有其他活跃的Warp可以切换，无法有效隐藏这部分巨大的内存延迟，因为几乎所有的Warp都因等待数据返回而陷入阻塞（Stall）状态。
优化思路：使用共享内存进行分块计算（Tiling）
#include<iostream> #include<vector> #include<random> #include<stdexcept> #include<cmath> #include<cuda_runtime.h> #define CUDA_CHECK(err) { \ cudaError_t err_ = (err); \ if (err_ != cudaSuccess) { \ std::cerr << "CUDA error at " << __FILE__ << ":" << __LINE__; \ std::cerr << ": " << cudaGetErrorString(err_) << std::endl; \ exit(EXIT_FAILURE); \ } \ } #define TILE_SIZE 32 __global__ void gemm_gpu(constfloat* A, constfloat* B, float* C, int M, int N, int K){ __shared__ float sA[TILE_SIZE][TILE_SIZE]; __shared__ float sB[TILE_SIZE][TILE_SIZE]; int row = blockIdx.y * TILE_SIZE + threadIdx.y; int col = blockIdx.x * TILE_SIZE + threadIdx.x; float sum = 0.0; for (int i = 0; i < (K+TILE_SIZE-1)/TILE_SIZE; ++i) { if (row <M && i*TILE_SIZE+threadIdx.x <K) sA[threadIdx.y][threadIdx.x]=A[row*K+i*TILE_SIZE+threadIdx.x]; else sA[threadIdx.y][threadIdx.x]=0.0; if ( i*TILE_SIZE+threadIdx.y <K && col <N) sB[threadIdx.y][threadIdx.x]=B[(i*TILE_SIZE+threadIdx.y)*N+col]; else sB[threadIdx.y][threadIdx.x]=0.0; __syncthreads(); for (int k =0; k< TILE_SIZE; ++k) { sum += sA[threadIdx.y][k]*sB[k][threadIdx.x]; } __syncthreads(); } if (row < M && col < N) { C[row * N + col] = sum; } } voidverify_result(constfloat* C_gpu, constfloat* C_cpu, int M, int N){ double max_error = 0.0; for (int i = 0; i < M * N; ++i) { max_error = fmax(max_error, fabs(C_gpu[i] - C_cpu[i])); } std::cout << "Max error between GPU and CPU results: " << max_error << std::endl; if (max_error > 1e-4) { std::cout << "Verification FAILED!" << std::endl; } else { std::cout << "Verification PASSED!" << std::endl; } } voidgemm_cpu(constfloat* A, constfloat* B, float* C, int M, int N, int K){ for (int r = 0; r < M; ++r) { for (int c = 0; c < N; ++c) { float sum = 0.0f; for (int i = 0; i < K; ++i) { sum += A[r * K + i] * B[i * N + c]; } C[r * N + c] = sum; } } } intmain(int argc, char** argv){ int M = 1024; int N = 1024; int K = 1024; std::cout << "Matrix dimensions: C(" << M << "x" << N << ") = A(" << M << "x" << K << ") * B(" << K << "x" << N << ")" << std::endl; std::vector<float> h_A(M * K); std::vector<float> h_B(K * N); std::vector<float> h_C_gpu(M * N, 0); std::vector<float> h_C_cpu(M * N, 0); std::random_device rd; std::mt19937 gen(rd()); std::uniform_real_distribution<> dis(0.0, 1.0); for (int i = 0; i < M * K; ++i) h_A[i] = dis(gen); for (int i = 0; i < K * N; ++i) h_B[i] = dis(gen); float *d_A, *d_B, *d_C; CUDA_CHECK(cudaMalloc(&d_A, M * K * sizeof(float))); CUDA_CHECK(cudaMalloc(&d_B, K * N * sizeof(float))); CUDA_CHECK(cudaMalloc(&d_C, M * N * sizeof(float))); CUDA_CHECK(cudaMemcpy(d_A, h_A.data(), M * K * sizeof(float), cudaMemcpyHostToDevice)); CUDA_CHECK(cudaMemcpy(d_B, h_B.data(), K * N * sizeof(float), cudaMemcpyHostToDevice)); dim3 block_dim(TILE_SIZE, TILE_SIZE); dim3 grid_dim((N + block_dim.x - 1) / block_dim.x, (M + block_dim.y - 1) / block_dim.y); std::cout << "\n--- Running GEMM ---" << std::endl; gemm_gpu<<<grid_dim, block_dim>>>(d_A, d_B, d_C, M, N, K); CUDA_CHECK(cudaDeviceSynchronize()); CUDA_CHECK(cudaMemcpy(h_C_gpu.data(), d_C, M * N * sizeof(float), cudaMemcpyDeviceToHost)); gemm_cpu(h_A.data(), h_B.data(), h_C_cpu.data(), M, N, K); verify_result(h_C_gpu.data(), h_C_cpu.data(), M, N); CUDA_CHECK(cudaFree(d_A)); CUDA_CHECK(cudaFree(d_B)); CUDA_CHECK(cudaFree(d_C)); return0; }
为了解决非合并访问与全局内存访问延迟高的问题，引入了基于共享内存的分块（Tiling）优化策略，其核心思想如图所示：
1. 任务划分：将输出矩阵C划分为若干个TILE_SIZExTILE_SIZE的数据块（Tile），并让一个线程块负责计算其中一个Tile。
2. 数据迭代加载与复用：如图要计算C的一个Tile，需要A矩阵的一个TILE_SIZE x K的行带和B矩阵的一个K x TILE_SIZE的列带。由于共享内存容量有限，无法一次性载入所有需要的数据，因此将K维度也进行分块，计算过程变为一个循环，每次迭代处理一个TILE_SIZExTILE_SIZE的A子矩阵和一个TILE_SIZExTILE_SIZE的B子矩阵。
代码时序图：
1. 加载A和B数据：从全局内存中将A和B的当前数据块加载到共享内存，这里的加载实现了内存合并访问，最大化内存吞吐。
2. 计算：在数据加载到共享内存后再并行地从共享内存中读取数据进行矩阵乘法计算，由于访问的是极低延迟的共享内存，计算速度非常快。
性能验证
全局内存事务次数对比
由前可知，非合并内存访问、全局内存访问都会增加全局内存事务次数，那么这里对比两次程序运行过程中执行的内存事务次数，确认次数量级的差异，这可以通过Nsight Compute办到。
==优化前== ---------------------------------------------- ----------- ------------- Metric Name Metric Unit Metric Value ---------------------------------------------- ----------- ------------- l1tex__t_bytes_pipe_lsu_mem_global_op_ld.sum byte 5,368,709,120 l1tex__t_bytes_pipe_lsu_mem_global_op_st.sum byte 4,194,304 l1tex__t_sectors_pipe_lsu_mem_global_op_ld.sum sector 167,772,160 l1tex__t_sectors_pipe_lsu_mem_global_op_st.sum sector 131,072 ---------------------------------------------- ----------- ------------- ==优化后== ---------------------------------------------- ----------- ------------ Metric Name Metric Unit Metric Value ---------------------------------------------- ----------- ------------ l1tex__t_bytes_pipe_lsu_mem_global_op_ld.sum byte 268,435,456 l1tex__t_bytes_pipe_lsu_mem_global_op_st.sum byte 4,194,304 l1tex__t_sectors_pipe_lsu_mem_global_op_ld.sum sector 8,388,608 l1tex__t_sectors_pipe_lsu_mem_global_op_st.sum sector 131,072 ---------------------------------------------- ----------- ------------
根据Nsight Compute统计的ld.global（从全局内存加载数据）指令和st.global（从全局内存写数据）指令上一共执行的事务次数，原来ld.global共计167772160次、st.global共计131072次（每次读/写内存事务处理32个字节），优化后ld.global共计8388608次、st.global共计131072次，从全局内存读取数据时执行的内存事务次数只有原来的5%，从全局内存读取的数据量从原来的5G降低为0.25G，说明更多的内存寻址offload到了更上层内存中。
程序执行耗时对比
程序耗时对比通过Nsight Compute也可以办到的，这里我用PTX插桩和Neutrino框架采集两个时间指标来验证性能改善情况，为后续章节抛砖引玉，这两个时间分别含义：
1. 线程在SM处理器上平均运行时间（Average Running Time）：涵盖指令发射+算术运算+内存访问整体时间
2. SM上平均调度空隙时间（Average Idle Time）：每个warp在等待调度执行时，因前序任务未及时释放资源而导致的平均等待时间（Warp stall、Warp发散或者Occupancy配置不当都有可能引发该时间变长）
具体采集实现细节后章节展开，这里只看结果，可以发现在经优化后的程序在两个时间上的耗时都缩短了，比较直观的可以看到优化对耗时整体提升是有效果的。
==优化前== Total Blocks: 1376, Total Tasks: 33024 Average Running Time: 190207 (cycles) Average Idle Time: 1196 (cycles) ==优化后== Total Blocks: 1024, Total Tasks: 32768 Average Running Time: 106999 (cycles) Average Idle Time: 474 (cycles)
持续优化
继续通过Nsight Compute进行性能分析，已经没有之前的非合并访问与全局内存访问频繁引发的性能问题了，新的提示是当前程序在写共享内存阶段存在bank conflict。
减少bank conflict
GPU硬件为了支持warp内32个线程对共享内存的并行访问，将共享内存组织为32个独立的存储体（banks）。一般理解上如果一个Warp内32个线程访问了同一个bank里不同地址，会引发bank conflict（如图所示），但是示例代码在一个Warp内32个线程的写入共享内存时应不存在bank conflict，这里提示存在bank conflict困扰了我挺久。
所以对bank conflict我有一个个人的不成熟理解，如果我把思考角度放到同一个block内的不同Warp之间，同指令所有线程并发访问时，是会存在不同行线程访问同一个bank的不同地址的情况，这或许是一种广义的bank资源争抢引发的bank conflict。
不管怎么样，因为每个bank在一个cycle内只能服务一次访问，因此这些冲突的访问就只能串行执行，带来了性能下降。
按照我的理解，对于示例中sB，它的每一列都访问了同一个bank号，引发了bank conflict：
这里我选择调整TILE_SIZE大小，从32减少至16，通过降低bank conflict的概率来进一步优化性能。
调整TILE_SIZE大小后bank conflict情况得到缓解，这里我同样编写了一个cuBLAS库代码(M,N,K=1024且fp32）进一步对比bank conflict，可见在将TILE_SIZE调整为16后的效果和使用cuBLAS实现的gemm在bank conflict层面量级差不多，这步优化暂时到此为止。
==TILE_SIZE 32== -------------------------------------------------------- ----------- ------------ Metric Name Metric Unit Metric Value -------------------------------------------------------- ----------- ------------ l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.avg 8,733.29 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.max 9,932 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.min 8,103 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.sum 628,797 -------------------------------------------------------- ----------- ------------ ==TILE_SIZE 16== -------------------------------------------------------- ----------- ------------ Metric Name Metric Unit Metric Value -------------------------------------------------------- ----------- ------------ l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.avg 2,309.79 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.max 2,745 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.min 1,774 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.sum 166,305 -------------------------------------------------------- ----------- ------------ ==cuBlAS== -------------------------------------------------------- ----------- ------------ Metric Name Metric Unit Metric Value -------------------------------------------------------- ----------- ------------ l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.avg 2,275.56 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.max 2,304 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.min 2,048 l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_st.sum 163,840 -------------------------------------------------------- ----------- ------------
小结
经过再次调整后的示例代码在使用的指令级别依然还有可持续优化空间，这里不继续展开，免得内容过于冗长。本部分内容主要是介绍CUDA基础，以及如何对Kernel进行性能分析，可以看到对Kernel性能分析上还是以Nsight Compute为主进行。不过本文并非讲Nsight Compute要如何分析，下面章节的内容，会开始介绍PTX的理论知识，如何进行插桩，以及如何对抓取结果进行可视化展示，介绍一些与Nsight Compute差异性的东西。
第二部分：从CUDA到PTX，看懂编译器生成的“中间语言”
1. PTX简介
不同代系GPU的原生机器码（SASS）不同，指令集存在差异，PTX是介于CUDA代码和底层硬件机器码之间的中间汇编语言，目的是提供公开、稳定、向前兼容的虚拟指令集架构，它可以在不同架构之间移植。
对于将在GPU上运行的程序，GPU编译器如nvcc会生成PTX，这段PTX汇编语言再通过一个名为ptxas的汇编程序生成可执行二进制文件（.cubin)，该二进制文件包含了为特定GPU架构编译号的SASS机器码，是GPU硬件直接解码和执行的二进制指令：
2. 看懂PTX
有几个途径可以看到对应CUDA代码的PTX汇编语言
1. 使用nvcc --keep -c 对CUDA C++文件编译时会保留.ptx文件
2. 使用cuobjdump -ptx将nvcc编译后的执行程序反汇编为PTX汇编语言
3. https://godbolt.org/此网站可以对CUDA C++代码在线提供PTX，高亮展示C++代码条目与汇编段落联系，初学者可参考：
这里依然以之前的示例程序初识PTX，下面是gemm.cu文件生成的PTX汇编内容，结合示例代码手工标注了关键汇编语言的含义，帮助看官们能看懂PTX的内容：
.version 8.7 .target sm_52 .address_size 64 --> 全局声明，标识这段代码的PTX版本，目标架构与地址模式 .visible .entry _Z8gemm_gpuPKfS0_Pfiii( --> 定义一个可被CPU调用的Kernel函数 .param .u64 _Z8gemm_gpuPKfS0_Pfiii_param_0, --> 参数声明 A .param .u64 _Z8gemm_gpuPKfS0_Pfiii_param_1, --> 参数声明 B .param .u64 _Z8gemm_gpuPKfS0_Pfiii_param_2, --> 参数声明 C .param .u32 _Z8gemm_gpuPKfS0_Pfiii_param_3, --> 参数声明 M .param .u32 _Z8gemm_gpuPKfS0_Pfiii_param_4, --> 参数声明 N .param .u32 _Z8gemm_gpuPKfS0_Pfiii_param_5 --> 参数声明 K ) { .reg .pred %p<9>; --> 寄存器声明9个bool寄存器 .reg .f32 %f<30>; --> 寄存器声明30个float寄存器 .reg .b32 %r<32>; --> 寄存器声明32个通用32位寄存器 .reg .b64 %rd<34>; --> 寄存器声明34个通过64位寄存器 ld.param.u64 %rd18, [_Z8gemm_gpuPKfS0_Pfiii_param_0]; --> 参数加载 A->%rd18 ld.param.u64 %rd19, [_Z8gemm_gpuPKfS0_Pfiii_param_1]; --> 参数加载 B->%rd19 ld.param.u64 %rd17, [_Z8gemm_gpuPKfS0_Pfiii_param_2]; --> 参数加载 C->%rd17 ld.param.u32 %r14, [_Z8gemm_gpuPKfS0_Pfiii_param_3]; --> 参数加载 M->%rd14 ld.param.u32 %r12, [_Z8gemm_gpuPKfS0_Pfiii_param_4]; --> 参数加载 N->%rd12 ld.param.u32 %r13, [_Z8gemm_gpuPKfS0_Pfiii_param_5]; --> 参数加载 K->%rd13 cvta.to.global.u64 %rd1, %rd19; --> 将通用地址转换为全局内存地址空间的64位指针得B基地址B[0][0] cvta.to.global.u64 %rd2, %rd18; mov.u32 %r15, %ntid.y; --> %ntid.y=blockDim.y,传入%r15 mov.u32 %r16, %ctaid.y; --> %ctaid.y=blockIdx.y,传入%r16 mov.u32 %r17, %tid.y; --> %tid.y=threadIdx.y,传入%r17 mad.lo.s32 %r1, %r16, %r15, %r17; --> %r16*%r15+%r17=row,计算结果传入%r1 mov.u32 %r18, %ntid.x; mov.u32 %r19, %ctaid.x; mov.u32 %r20, %tid.x; mad.lo.s32 %r2, %r19, %r18, %r20; setp.ge.s32 %p1, %r1, %r14; --> 此处用的是%r1是否>=%r14来判断row<M是否成立，判断结果放入%p1 setp.ge.s32 %p2, %r2, %r12; or.pred %p3, %p1, %p2; --> 编译器优化，判断(row>=M||col>=N) @%p3 bra $L__BB0_9; --> 判断%p3结果是否True，如果真跳转到$L__BB0_9(就是ret结束） setp.lt.s32 %p4, %r13, 1; --> 编译器优化，用%r13是否<1判断K是否=0,判断K是否=0 mov.f32 %f29, 0f00000000; --> 根据代码上下文可知，%f29=sum，赋值0.0f @%p4 bra $L__BB0_8; --> 如果K=0goto $L__BB0_8，不会进入循环体 add.s32 %r22, %r13, -1; --> K-1的结果传入%r22 and.b32 %r31, %r13, 3; --> 按位与计算 K&3，等价于K%4（4路展开是由编译器控制，表示一次循环控制每次迭代处理4次原循环的工作，代码里可以通过#pragma unroll独立指定 setp.lt.u32 %p5, %r22, 3; --> 判断%r22是否<3,等价于判断K是否小于4 mov.f32 %f29, 0f00000000; mov.u32 %r30, 0; --> %r30是i，这里就是给i赋值0 @%p5 bra $L__BB0_5; --> 如果K<4为True，goto $L__BB0_5 sub.s32 %r29, %r13, %r31; --> （4路展开初始化部分）K-（K%4）传入%r29，这里是判断K是否可以被4整除，用于后续不能被4整除的部分goto $L__BB0_5。 计算得到的%r29是能被4整除部分，用于后续控制4路循环循环次数 mul.lo.s32 %r24, %r13, %r1; --> K*row 计算 A[row][i] 的行偏移，传入%r24 mul.wide.s32 %rd3, %r24, 4; --> 计算A[row][0]的字节偏移 mul.wide.s32 %rd20, %r2, 4; --> 以行首为基准，计算第col列的字节偏移，即计算B[0][col]的字节偏移 add.s64 %rd30, %rd1, %rd20; --> 基地址+字节偏移=B[0][col]的起始地址 mul.wide.s32 %rd5, %r12, 4; --> N*4=B矩阵的每行字节数（stride） mov.f32 %f29, 0f00000000; mov.u32 %r30, 0; mov.u64 %rd31, %rd2; --> A基地址传入%rd31 $L__BB0_4: --> K > 4, 这里是示例代码配置下实际会执行的循环体部分 add.s64 %rd21, %rd31, %rd3; --> 基地址+字节偏移=A[row][0]起始地址 ld.global.f32 %f12, [%rd30]; --> 从全局内存地址%rd30处读取一个float，放入寄存器%f12 ld.global.f32 %f13, [%rd21]; fma.rn.f32 %f14, %f13, %f12, %f29; --> %f13*%f12+%f29，对于首次循环就是计算sum=A[row][0]*B[0][col]+sum add.s64 %rd22, %rd30, %rd5; --> B[1][col] ld.global.f32 %f15, [%rd22]; ld.global.f32 %f16, [%rd21+4]; --> A[row][1] fma.rn.f32 %f17, %f16, %f15, %f14; add.s64 %rd23, %rd22, %rd5; --> B[2][col] ld.global.f32 %f18, [%rd23]; ld.global.f32 %f19, [%rd21+8]; --> A[row][2] fma.rn.f32 %f20, %f19, %f18, %f17; add.s64 %rd24, %rd23, %rd5; add.s64 %rd30, %rd24, %rd5; ld.global.f32 %f21, [%rd24]; ld.global.f32 %f22, [%rd21+12]; fma.rn.f32 %f29, %f22, %f21, %f20; add.s32 %r30, %r30, 4; --> 4路展开所以下一循环从i+4开始 add.s64 %rd31, %rd31, 16; add.s32 %r29, %r29, -4; --> 循环次数减4 setp.ne.s32 %p6, %r29, 0; --> 判断%r29是否减到0了 @%p6 bra $L__BB0_4; --> 如果不等于0，goto $L__BB0_4 开始为止继续计算 $L__BB0_5: --> K<4 ，普通单路循环，处理掉不够4次的情况 setp.eq.s32 %p7, %r31, 0; @%p7 bra $L__BB0_8; mad.lo.s32 %r25, %r30, %r12, %r2; mul.wide.s32 %rd25, %r25, 4; add.s64 %rd33, %rd1, %rd25; mul.wide.s32 %rd11, %r12, 4; mad.lo.s32 %r26, %r13, %r1, %r30; mul.wide.s32 %rd26, %r26, 4; add.s64 %rd32, %rd2, %rd26; $L__BB0_7: .pragma "nounroll"; ld.global.f32 %f23, [%rd33]; ld.global.f32 %f24, [%rd32]; fma.rn.f32 %f29, %f24, %f23, %f29; add.s64 %rd33, %rd33, %rd11; add.s64 %rd32, %rd32, 4; add.s32 %r31, %r31, -1; setp.ne.s32 %p8, %r31, 0; @%p8 bra $L__BB0_7; $L__BB0_8: --> 对应到是 C[row * N + col] = sum mad.lo.s32 %r27, %r1, %r12, %r2; --> 对应 row * N + col cvta.to.global.u64 %rd27, %rd17; mul.wide.s32 %rd28, %r27, 4; --> %r27*4字节，传入%rd28，这是在计算C元素的内存地址 add.s64 %rd29, %rd27, %rd28; st.global.f32 [%rd29], %f29; --> 写入全局内存 $L__BB0_9: ret; }
从PTX的内容也可以看出，原示例里多是对global memory（ld.global和st.global）的操作，而且Line 57 这里%rd5的值特别大（1024*4字节），在读取B元素的时候每次需要跨4096个字节。
3. 小结
这部分主要介绍什么是PTX，以及如何看懂一个简单的PTX汇编内容，有助于理解后续章节基于PTX插桩和分析的内容。
第三部分：深入 PTX，构建可编程的GPU性能探针
1. 为什么选择PTX
从下图中可以看出来，无论是CUDA C++预编译，还是市面上主流的大模型框架（PyTorch、TensorFlow）等基于JIT即时编译的模式，它们在针对GPU时，代码生成层都是相同的PTX代码，这是连接上层软件生态和下层多代GPU硬件的纽带，基于PTX可以构建通用的分析工具：
2. GPU 上的“eBPF 式”探测模型
我借助的是Neutrino这个可编程探针框架来实现我的需求，它利用在PTX层插桩来实现对指令级的细粒度性能观测，支持从Warp/Thread两个级别进行各项自定义分析任务，观测时对实际性能开销较低。
框架介绍
从宏观来看，它和eBPF有相似之处，如下图，性能“探针”（Probe）是由三个核心部分组合而成的：要执行的代码（Snippet）、执行的时机（Tracepoint）和存放结果的数据结构（Map）：
- Snippet是一小段PTX汇编代码，定义这个探针要做什么，功能类似于eBPF里某个kprobe具体执行的代码体，只是代码内容从C语言换做PTX汇编语言；
- Tracepoint定义Snippet代码被触发执行的时机；
- Map定义了用于存储探测结果的内存布局与数据结构，功能也于eBPF的map相似；
以下这个图，详细解释了Neutrino这个工具探针执行隔离和数据存储机制，使用蓝色代表原始Kernel代码，紫色代码探针代码：
- 在指令层面（图A区），探针指令被“见缝插针”地插入到原始指令流中，它们在同一个时间线上执行，但逻辑上是独立的。
- 在寄存器层面（图A区），原始代码和探针代码使用的寄存器做了隔离，两套寄存器完全分离确保不会意外地覆盖或破坏原始程序正在使用的重要数据，保证探针的无侵入性。
- 在内存层面（图B区），Neutrino框架在全局内存中预先分配的一块结构化缓冲区，这是一个二维网格，每个线程都有自己在Probe Map中唯一、专属的存储位置，每个线程都只往自己的格子里写数据，不会发生写入冲突。
下图是整个框架工作流程，解释了框架是如何通过Hooking GPU驱动、即时编译探针和动态插桩来实现对任意GPU程序的透明化性能分析的。可以把整个流程看作是三大利益相关方的互动：用户（绿色部分）、被分析的程序（黑色/浅蓝色部分）和Neutrino的核心引擎（蓝色和紫色部分）：
- 用户通过
neutrino -p dmat python main.py方式启动探针和被分析程序，-p dmat是预定义的探针名称，python main.py是用户想要分析的目标程序 - Neutrino在运行目标程序前（绿色框），从Py DSL Probes--通过Python编写的DSL加载自定义探针内容，编译器会把DSL探针定义转换为以TOML包装的GPU汇编语言（ASM）探针，进行一系列的Verifier安全验证器的检查后（比如语义检查、是否覆写寄存器合法地址），探针代码会被设置到一个环境变量以便后续组件访问
- 准备工作完成，工具使用fork和exec启动用户指定的目标程序，它通过LD_PRELOAD技术在目标程序启动时，强制加载一个自定义的动态链接库，这个库就是Hooked Driver
- （蓝色框）这个Hooked Driver取代了系统原始的GPU驱动的部分功能，会拦截所有目标程序发出的GPU相关的API调用，调用Probe Engine把拦截的code交给它处理
- （紫色框）Probe Engine先将交给它处理的Kernel二进制代码反编译为汇编指令，然后根据之前从环境变量读取的探针定义，将探针的汇编指令插入到原始汇编代码的正确位置，最后再把插桩后的新汇编代码重新组装成一个可执行的二进制probed code
- （返回蓝色框）Hooked Driver收到Probe Engine返回的probed code，将插桩后的代码、原始输入数据的地址和Probe map的地址一起传递给GPU后，真正调用GPU驱动API在GPU上执行这个被修改过的Kernel。待GPU执行完毕，Hooked Driver将Probe map中收集到的数据从GPU内存转储到磁盘上用于离线分析
步骤拆解
编写Probe
实际一个Probe探针代码如下所示，Neutrino提供了DSL作为探测器的高级接口，类似eBPF的bpftrace，不过DSL有一定功能局限性，它不是唯一的选择，也是可以直接编写汇编语言来实现高级用途，下一个小章节会通过示例展示直接通过编写汇编来实现的自定义探针。
我先解析一下这个Python形式的探针实现的含义（本文第一部分统计程序执行耗时，就是用的这个探针）。
from neutrino import probe, Map import neutrino.language as nl CALLBACK = "block_sched_callback.py"# for trace analysis # declare maps for persistence @Map(level="warp", type="array", size=16, cap=1) class block_sched: start: nl.u64 elapsed: nl.u32 cuid: nl.u32 # declare probe registers shared across probes start: nl.u64 = 0# starting clock elapsed: nl.u64 = 0# elapsed time, initialized to 0 # define probes with decorator @probe(pos="kernel", level="warp", before=True) def thread_start(): start = nl.clock() @probe(pos="kernel", level="warp") def thread_end(): elapsed = nl.clock() - start block_sched.save(start, elapsed, nl.cuid())
- Line1-2：导入Neutrino的相关依赖包
- Line4：指定一个回调脚本，这是用来对收集到的数据进行按需分析的Python文件
- Line7-11：定义Map数据存储结构，@Map(...)装饰器定义Probe Map的元数据，其中level="warp"，表示数据存储粒度是每个Warp，type="array", size=16, cap=1 是为每个Warp创建一个容量为1的数组（仅采集Warp里Leader Thread），包含16字节（8个字节start、4字节elapsed、4字节cuid）。class block_sched则是定义最终存储在磁盘的数据内容，每个记录都会包含start、elapsed、cuid
- Line14-15：初始化probe regs
- Line18-20：@probe()装饰器定义Tracepoint，其中pos="kernel", level="warp", before=True表示跟踪点是整个Kernel范围，定义的Snippet要插桩到kernel开始执行前（不加before默认是after）。def函数体定义具体的Snippet代码，start = nl.clock()会被编译为一条汇编指令：mov.u64 %NR0, %clock64，将kernel执行前的时钟值存入Neutrino额外声明的%NR0寄存器。
- Linux22-25：定义第二个Tracepoint，目的在kernel执行后插桩，elapsed = nl.clock() - start会被编译为：mov.u64 %NR2, %clock64; sub.u64 %NR1, %NR2, %NR0; 即在计算时间差。block_sched.save()是Neutrino提供的辅助函数，表示把这些值保存到名为block_sched到Map中，同样save指令也是编译成汇编指令执行，比如通过@%leader st.global.v2.u32 [%map_block_sched1], { %block_sched_u32_1, %block_sched_u32_2 } 存入全局内存。
执行Probe
当探针hook到目标程序执行后，会生成一个trace目录：
[info] trace in ./trace/Aug14_112435_2682
目录包含如下文件：
./trace/Aug14_112435_2682 ├── event.log --> 工具的事件日志 ├── kernel --> 目标程序里每个Kernel会生成一个子目录 │ ├── 0_d02df43fa8dfe9acebcd959bee6f98c2404e40f9 │ │ ├── kernel.info --> 记录Kernel名称、形状 │ │ ├── original.bin --> 原Kernel的cubin文件 │ │ ├── original.ptx --> 原Kernel的ptx汇编文件 │ │ ├── probed.bin --> 插桩后生成的cubin文件 │ │ ├── probed.ptx --> 插桩后的ptx汇编文件 │ │ ├── process.log --> ptxas对插桩后的ptx生成cubin的日志记录 │ │ ├── pruned.bin │ │ └── pruned.ptx --> 经过裁剪的仅包含与探针相关逻辑的PTX代码片段 │ └── 1_f3396b6dedff3f84b1c26dccb7524fb932dd6c14 │ ├── kernel.info │ ├── original.bin │ ├── original.ptx │ ├── probed.bin │ ├── probed.ptx │ ├── process.log │ ├── pruned.bin │ └── pruned.ptx ├── probe.toml --> toml格式的探针文件 ├── read.py --> Neutrino预置的用来解析result目录中二进制数据的脚本 └── result ├── 0.297500.bin --> 通过Map收集到的所有原始二进制数据,每个Kernel一个 └── 5.111352.bin
插桩后生成的ptx不再逐条讲解，各位看官可以看图对比差异：
数据分析
然后就到了最后，对最终抓取的二进制结果提炼有意义的分析数据，这步分为两部分：
- 解析二进制数据文件（例如0.297500.bin），将其还原为结构化，人类可读的记录；
- 基于解析的数据再进行数据分析；
解析数据
阅读Neutrino源码可以了解到.bin存储时遵循特定的自定义结构，由三部分组成（如图），在自定义探针Map数据前还有两个预置传入的类，读取时需要按照该顺序读取数据，避免数据错乱。
使用Python进行二进制文件解析，代码如下面translate.py，即可得到捕获的探测数据，这看到的数据是每个block里每个warp中leader thread的开始时间、持续时间、所在SM处理器编号。
importstruct # i表示int32，4字节 fmt_header ='iiiiiiii' header_size = struct.calcsize(fmt_header) fmt_section = 'iii' section_size = struct.calcsize(fmt_section) # q表示int64，8字节 fmt_data = 'qii' data_size = struct.calcsize(fmt_data) events = [] with open("0.297500.bin", "rb") as f: #读取&丢弃 header_data = f.read(header_size) #读取&丢弃 section_data = f.read(section_size) offset = struct.unpack(fmt_section, section_data)[2] #跳转到自定义Map数据起始位置 f.seek(offset) while chunk := f.read(data_size): if len(chunk) == data_size: start, elapsed, cuid = struct.unpack(fmt_data, chunk) events.append((start, elapsed, cuid)) print(f"start={start}, elapsed={elapsed}, cuid={cuid}")
Neutrino提供了预置的用于解析二进制的脚本read.py（在生成的trace目录里），读取原理和我自己写的translate.py一致，比如block_sched这个探针的read.py内容如下：
处理数据
对数据的处理是根据我们实际的目的来进行的，最终数据展示可以是多种多样。
第一部分里讨论的程序执行耗时，两个时间的处理完整实现代码如下：
# Neutrino Auto-Generated Code for Trace Reading importstruct fromtypingimportNamedTuple, List, Tuple, Dict fromneutrinoimportTraceHeader, TraceSection classblock_sched(NamedTuple): start: int elapsed: int cuid: int def parse(path: str) -> Tuple[TraceHeader, List[TraceSection], Dict[str, List[List[NamedTuple]]]]: with open(path, "rb") as f: header: TraceHeader = TraceHeader(*struct.unpack("iiiiiiii", f.read(32))) sections: List[TraceSection] = [] for _ in range(header.numProbes): sections.append(TraceSection(*struct.unpack("IIQ", f.read(16)))) gridSize = header.gridDimX * header.gridDimY * header.gridDimZ blockSize = header.blockDimX * header.blockDimY * header.blockDimZ records: Dict[str, List[List[NamedTuple]]] = dict() # Read block_sched records["block_sched"] = [] f.seek(sections[0].offset) for i in range(gridSize): records["block_sched"].append([]) for j in range(blockSize // sections[0].warpDiv): records["block_sched"][-1].append([]) for k in range(sections[0].size // 16): records["block_sched"][i][j].append(block_sched(*struct.unpack("qII", f.read(16)))) return header, sections, records # END of Neutrino Auto-Generated Code for Trace Reading import sys import numpy as np from typing import NamedTuple, List header, sections, records_map = parse(sys.argv[1]) # filled by path to trace records = records_map["block_sched"] all_tasks = [task for block in records for warp in block for task in warp] all_tasks.sort(key=lambda task: task.start) unique_sms = sorted(list(set(task.cuid for task in all_tasks))) sm_to_idx = {sm_id: i for i, sm_id in enumerate(unique_sms)} num_sms = len(unique_sms) sm_timelines = [[] for _ in range(num_sms)] sched_times = [0.0] * num_sms work_times = [0.0] * num_sms for task in all_tasks: sm_idx = sm_to_idx[task.cuid] timeline = sm_timelines[sm_idx] #按照SM处理器粒度，统计线程耗时总和 work_times[sm_idx] += task.elapsed ifnot timeline: timeline.append(task) continue best_slot = None earliest_end_time = float('inf') for old_task in timeline: if old_task.start + old_task.elapsed <= task.start: if old_task.start + old_task.elapsed < earliest_end_time: earliest_end_time = old_task.start + old_task.elapsed best_slot = old_task if best_slot is not None: sched_delay = task.start - (best_slot.start + best_slot.elapsed) #按照SM处理器粒度，统计Warp空闲时间 sched_times[sm_idx] += sched_delay timeline.remove(best_slot) timeline.append(task) else: timeline.append(task) avg_work_time = np.array(work_times).sum() / len(all_tasks) if all_tasks else0 avg_idle_gap = np.array(sched_times).sum() / len(all_tasks) if all_tasks else0 print(f"Total Blocks: {len(records)}, Total Tasks: {len(all_tasks)}") print(f"Average Running Time: {int(avg_work_time)} (cycles)") print(f"Average Idle Time: {int(avg_idle_gap)} (cycles)")
题外话
Probe探测纬度分Warp和Thread两个，Warp主要用于耗时分析，Thread主要用于值分析。
前面使用的Probe探针插桩位置在Kernel，也支持在具体的ptx指令的前后插桩，比如下图，分别在读取全局内存和写入全局内存的指令后面进行插桩，统计每个线程的全局内存总流量（不过这个统计数据不准确，ld.global会经过L1 Cache，实际并不会全都访问全局内存），这结果更加倾向于是用于侧面参考。
写在本文里是做个插桩位置在指令层面的效果展示。
3. 示例：捕获每个Warp执行耗时
回到引言，比如某个Kernel执行时第 128 个warpid到底花了多少时间，how to do it？这里目标程序采用示例gemm.cu。
探针的逻辑是根据ptx特殊寄存器获取threadIdx、blockIdx、blockDim等信息，计算出全局范围内的warpIdx，统计每个warp的leader thread用时，这样就获取到每个warp耗时，完整的探针代码和数据分析代码如下，通过neutrino -p warp_duration.toml ./gemm > warp_duration.log执行。
regs = 25 CALLBACK = "warp_duration_callback.py" [probe.warp_duration_trace] level = "warp" pos = "kernel" before = """ .reg .b32 %tmp<21>; .reg .pred %lead; mov.u32 %tmp2, %laneid; setp.eq.u32 %lead, %tmp2, 0; @%lead mov.u32 %tmp3, %nwarpid; @%lead mov.u32 %tmp4, %tid.x; // threadIdx.x @%lead mov.u32 %tmp5, %tid.y; // threadIdx.y @%lead mov.u32 %tmp6, %tid.z; // threadIdx.z @%lead mov.u32 %tmp7, %ntid.x; // blockDim.x @%lead mov.u32 %tmp8, %ntid.y; // blockDim.y @%lead mov.u32 %tmp18, %ntid.z; // blockDim.z @%lead mov.u32 %tmp9, %ctaid.x; // blockIdx.x @%lead mov.u32 %tmp10, %ctaid.y; // blockIdx.y @%lead mov.u32 %tmp11, %ctaid.z; // blockIdx.z @%lead mov.u32 %tmp12, %nctaid.x; // gridDim.x @%lead mov.u32 %tmp13, %nctaid.y; // gridDim.y @%lead mad.lo.s32 %tmp14, %tmp8, %tmp6, %tmp5; // blockDim.y * threadIdx.z + threadIdx.y @%lead mad.lo.s32 %tmp15, %tmp14, %tmp7, %tmp4; // thread_idx in per block = (blockDim.y * threadIdx.z + threadIdx.y) * blockDim.x + threadIdx.x @%lead div.s32 %tmp15, %tmp15, 32; // get persistent warpid instead of dynamic %warpid @%lead mad.lo.s32 %tmp16, %tmp13, %tmp11, %tmp10; // gridDim.y * blockIdx.z + blockIdx.y @%lead mad.lo.s32 %tmp17, %tmp16, %tmp12, %tmp9; // block_idx = (gridDim.y * blockIdx.z + blockIdx.y) * gridDim.x + blockIdx.x @%lead mul.lo.s32 %tmp19, %tmp7, %tmp8; @%lead mul.lo.s32 %tmp20, %tmp19, %tmp18; // thread_num=blockDim.x*blockDim.y*blockDim.z @%lead div.s32 %tmp20, %tmp20, 32; // warp_num per block = thread_num/32 @%lead mad.lo.s32 %tmp1, %tmp17, %tmp20, %tmp15; // buf_idx = block_idx * warp_num + warpIdx @%lead mov.u64 %NR0, %clock64;""" after = """ mov.u64 %NR2, %clock64; sub.u64 %NR1, %NR2, %NR0; cvt.u64.u32 %NR3, %tmp1; SAVE [ warp_duration ] { %NR0, %NR1, %NR3 }""" [map.warp_duration] level = "warp" type = "array" size = 32 cap = 1 [map.warp_duration.regs] start = [ "u64", "None",] elapsed = [ "u64", "None",] warpid = [ "u64", "None",]
importstruct fromtypingimportNamedTuple fromneutrinoimportTraceHeader, TraceSection classwarp_duration(NamedTuple): start: int elapsed: int warpid: int def parse(path: str) -> tuple[TraceHeader, list[TraceSection], dict[str, list[list[NamedTuple]]]]: with open(path, "rb") as f: header: TraceHeader = TraceHeader(*struct.unpack("iiiiiiii", f.read(32))) sections: list[TraceSection] = [] for _ in range(header.numProbes): sections.append(TraceSection(*struct.unpack("IIQ", f.read(16)))) gridSize = header.gridDimX * header.gridDimY * header.gridDimZ blockSize = header.blockDimX * header.blockDimY * header.blockDimZ warp_num = blockSize*gridSize // 32 records: dict[str, list[list[NamedTuple]]] = dict() records["warp_duration"] = [] warp_duration_section = sections[0] f.seek(warp_duration_section.offset) per_warprecord_size = warp_duration_section.size total_buffer_size = warp_num * per_warprecord_size num_entries = warp_num print(f"正在读取records, 每个记录大小={per_warprecord_size}, 总记录数={num_entries}") raw_data = f.read(total_buffer_size) for i in range(num_entries): offset = i * per_warprecord_size if offset + per_warprecord_size > len(raw_data): print(f"警告: 文件数据不完整，提前结束于第 {i} 条记录。") break start, elapsed, warpid = struct.unpack_from("QQQ", raw_data, offset) records["warp_duration"].append(warp_duration(start=start, elapsed=elapsed, warpid=warpid)) from collections import defaultdict grouped = defaultdict(list) for record in records["warp_duration"]: if record.start > 0or record.elapsed > 0: grouped[record.warpid].append(record) # 异常筛选 zero_id_records = [r for r in records["warp_duration"] if r.start == 0and r.elapsed == 0] sorted_warpids = sorted(grouped.keys()) print("\n--- Warp Execution Duration Summary ---") print(f"{'Warp ID':<10} {'Elapsed Cycles (avg)':<20} {'Count':<8} {'Start (first)'}") for wid in sorted_warpids: entries = grouped[wid] avg_elapsed = sum(e.elapsed for e in entries) / len(entries) valid_starts = [e.start for e in entries] first_start = min(valid_starts) if valid_starts else0 print(f"{wid:<10} {avg_elapsed:<20.2f} {len(entries):<8} {first_start}") if zero_id_records: print(f"\n注意: 发现了 {len(zero_id_records)} 条无效记录 (start, elapsed均为0)。") print("可能这些Warps的计时数据没有被成功写回。") return header, sections, records import sys import numpy as np from typing import NamedTuple, List header, sections, records_map = parse(sys.argv[1]) # filled by path to trace records = records_map["warp_duration"] #print(f"check point 1: original data: {records}")
采样结果通过tail -n +2 warp_duration.log | sort -k4 -n，按照warp开始执行时间从早到晚排序，可以发现在Kernel执行时，不是按照顺序从序号0开始调度Warp，并且由于L2/DRAM资源争抢/拥塞，呈现当越多的Warp启动，Warp耗时变长的现象，如图：
Warp ID Elapsed Cycles(avg) Count Start(first) 261 138413.00 1 107148575 262 138400.00 1 107148577 263 138408.00 1 107148579 260 138388.00 1 107148588 1125 137505.00 1 107148590 1126 137552.00 1 107148592 1127 137580.00 1 107148594 1988 138808.00 1 107148594 1989 138846.00 1 107148596 1990 138824.00 1 107148598 1991 138832.00 1 107148600 1124 137521.00 1 107148603 ..... ..... 8629 196405.00 1 107451009 7929 195593.00 1 107451010 8212 194960.00 1 107451011 8352 195822.00 1 107451012 7829 195584.00 1 107451013 7925 195658.00 1 107451014 8208 195607.00 1 107451015 8215 195331.00 1 107451022 8533 195631.00 1 107451022 8632 196351.00 1 107451033 .....
对比gemm_v2.cu探测数据，如下图，每个Warp耗时明显缩短：
4. 小结
通常情况下，在GPU性能调优实践中，Nsight Compute已经足够应对大多数场景。本文的PTX插桩方式更多情况下，是提供给CUDA开发者“内视”Kernel的技术，在更深层次的Kernel分析时，基于PTX插桩的方法就展现出其独特价值。
我作为一名GPU性能优化初学者，对底层的探索还很粗浅，但掌握了这种深入到指令层面的分析方法，无疑让我们在面对复杂的GPU性能“疑难杂症”时，多了一份底气和一把利器。基于PTX插桩赋予了我提出并解答自己问题的能力。
深入到PTX指令层，能够帮助我们最直观地看到指令的交互模式，正是我们从“知道性能差”到“理解为何差”的关键一步。
学无止境，精益求精，这趟深入Kernel的旅程，绝对物有所值。
参考
https://docs.nvidia.com/cuda/cuda-c-programming-guide/#
https://docs.nvidia.com/cuda/cuda-runtime-api/index.html
https://docs.nvidia.com/cuda/parallel-thread-execution/#
https://arxiv.org/pdf/2208.11174
https://siboehm.com/articles/22/CUDA-MMM
https://feldmann.nyc/blog/smem-microbenchmarks
https://www.usenix.org/system/files/osdi25-huang-songlin.pdf
https://open-neutrino.github.io/docs
https://github.com/open-neutrino/neutrino/blob/main/README.md
来源 | 阿里云开发者公众号
作者 | 晨染
