1. 從一張顯卡的“脾氣”說(shuō)起為什么你的CUDA程序跑不快很多人第一次接觸CUDA編程心態(tài)是這樣的CPU上跑得慢那就扔到GPU上幾千個(gè)核心一起算怎么著也得快個(gè)幾十倍吧。結(jié)果代碼寫完一跑發(fā)現(xiàn)比CPU還慢或者只快了兩三倍離預(yù)期差得遠(yuǎn)。這時(shí)候就開(kāi)始懷疑人生——是不是顯卡不行是不是驅(qū)動(dòng)沒(méi)裝好其實(shí)問(wèn)題往往不在硬件而在于你沒(méi)有理解GPU的“脾氣”。GPU和CPU的設(shè)計(jì)哲學(xué)完全不同。CPU像是一個(gè)博士什么復(fù)雜邏輯都能處理但一次只能干幾件事GPU像是一萬(wàn)個(gè)中學(xué)生每個(gè)人只會(huì)做簡(jiǎn)單算術(shù)但一萬(wàn)人同時(shí)動(dòng)手吞吐量就上來(lái)了。如果你讓這一萬(wàn)個(gè)中學(xué)生去做博士的活兒——比如復(fù)雜的條件分支、頻繁的內(nèi)存跳轉(zhuǎn)——他們反而會(huì)互相拖后腿。這篇文章要聊的就是怎么順著GPU的架構(gòu)特點(diǎn)去寫CUDA代碼讓它真正跑出該有的性能。我會(huì)從GPU的硬件架構(gòu)講起把SM、warp、shared memory這些概念用大白話拆開(kāi)然后落到具體的優(yōu)化手段上怎么合并內(nèi)存訪問(wèn)、怎么用shared memory做tiling、怎么避免bank conflict、怎么選block和grid的尺寸、怎么用nsight compute去定位瓶頸。每一部分都會(huì)給出可復(fù)現(xiàn)的代碼示例和實(shí)測(cè)數(shù)據(jù)不是紙上談兵。適合誰(shuí)看如果你寫過(guò)CUDA但性能不理想或者正準(zhǔn)備把某個(gè)計(jì)算密集型任務(wù)搬到GPU上又或者你只是想知道“為什么我的矩陣轉(zhuǎn)置這么慢”那這篇內(nèi)容應(yīng)該能幫到你。不需要你是并行計(jì)算專家但至少要能看懂C語(yǔ)言和基本的線性代數(shù)。2. GPU架構(gòu)到底長(zhǎng)什么樣把硬件拆開(kāi)來(lái)看2.1 從SM到warpGPU的“車間”和“班組”先建立一個(gè)直觀的模型。一塊GPU芯片上有很多個(gè)SMStreaming Multiprocessor流多處理器你可以把每個(gè)SM想象成一個(gè)車間。一個(gè)車間里有若干個(gè)warp scheduler調(diào)度器每個(gè)調(diào)度器管著一組warp。warp是什么32個(gè)線程捆在一起像一個(gè)班組必須同時(shí)執(zhí)行同一條指令。這就是所謂的SIMTSingle Instruction, Multiple Threads模型。為什么是32這是硬件設(shè)計(jì)的選擇不是隨便定的。NVIDIA從很早的架構(gòu)開(kāi)始就采用32作為warp size一直延續(xù)到現(xiàn)在。32個(gè)線程共享一個(gè)指令發(fā)射端口如果這32個(gè)線程都走同一條路徑那效率最高如果遇到if-else分支一部分線程走if一部分走else那就得串行執(zhí)行兩條路徑這叫warp divergencewarp分歧是性能殺手之一。每個(gè)SM里還有寄存器文件、shared memory/L1 cache、以及各種執(zhí)行單元FP32、INT32、Tensor Core等。寄存器是每個(gè)線程私有的速度最快但數(shù)量有限。shared memory是block內(nèi)共享的速度僅次于寄存器但需要手動(dòng)管理。global memory就是顯存容量大但延遲高通常有400-800個(gè)時(shí)鐘周期的延遲。理解這個(gè)層級(jí)關(guān)系很重要寄存器 shared memory L2 cache global memory。優(yōu)化的核心思路就是盡量讓數(shù)據(jù)在快的層級(jí)里被反復(fù)使用減少對(duì)慢層級(jí)的訪問(wèn)。2.2 占用率不是越高越好理解occupancy的真實(shí)含義Occupancy占用率是指一個(gè)SM上活躍warp數(shù)量與最大支持warp數(shù)量的比值。很多人以為occupancy越高越好拼命調(diào)block size想把occupancy拉滿。但實(shí)際上occupancy只是手段不是目的。高occupancy的好處是當(dāng)一個(gè)warp在等內(nèi)存的時(shí)候其他warp可以頂上隱藏延遲。但如果你每個(gè)線程用了大量寄存器occupancy自然就上不去。這時(shí)候如果強(qiáng)行降低寄存器用量可能會(huì)導(dǎo)致register spilling寄存器溢出到local memory反而更慢。我實(shí)測(cè)過(guò)一個(gè)向量加法的kernelblock size從128調(diào)到1024occupancy從50%拉到100%但性能只提升了不到5%。因?yàn)橄蛄考臃ū旧砭褪莔emory-bound瓶頸在顯存帶寬不在計(jì)算單元。相反在一個(gè)矩陣乘法的kernel里通過(guò)shared memory tiling把occupancy從25%提到50%性能翻了將近一倍因?yàn)闇p少了global memory訪問(wèn)。所以正確的做法是先用nsight compute看一下kernel是compute-bound還是memory-bound再?zèng)Q定優(yōu)化方向。不要盲目追求occupancy數(shù)字。2.3 內(nèi)存層級(jí)與訪問(wèn)延遲為什么合并訪問(wèn)如此關(guān)鍵Global memory的訪問(wèn)是以32字節(jié)為單位一個(gè)sector的。當(dāng)一個(gè)warp的32個(gè)線程訪問(wèn)連續(xù)的內(nèi)存地址時(shí)硬件可以把這些訪問(wèn)合并成少數(shù)幾個(gè)transaction。比如每個(gè)線程訪問(wèn)4字節(jié)float32個(gè)線程就是128字節(jié)正好是4個(gè)sector一次搞定。但如果線程訪問(wèn)的地址是跳躍的比如stride為2那128字節(jié)的數(shù)據(jù)只用了64字節(jié)浪費(fèi)了一半帶寬。這就是所謂的coalesced access合并訪問(wèn)。它是CUDA優(yōu)化里最基礎(chǔ)也最重要的一條規(guī)則。很多人寫kernel的時(shí)候習(xí)慣性地按照CPU的思維去索引數(shù)組結(jié)果導(dǎo)致stride訪問(wèn)帶寬利用率只有20%-30%。舉個(gè)例子矩陣按行存儲(chǔ)如果warp里的線程按列訪問(wèn)那每個(gè)線程的地址間隔就是矩陣的寬度完全無(wú)法合并。解決辦法要么是轉(zhuǎn)置矩陣要么用shared memory做中轉(zhuǎn)。后面會(huì)詳細(xì)講。3. CUDA編程優(yōu)化的核心手段從內(nèi)存到指令3.1 合并訪問(wèn)與向量化加載讓帶寬跑滿先看一段“反面教材”__global__ void badCopy(float* out, float* in, int n) { int idx threadIdx.x blockIdx.x * blockDim.x; if (idx n) { out[idx] in[idx]; } }這段代碼其實(shí)是合并訪問(wèn)的因?yàn)橄噜従€程訪問(wèn)相鄰地址。但如果我們把索引改成idx (threadIdx.x blockIdx.x * blockDim.x) * 2那就變成了stride-2訪問(wèn)帶寬直接減半。更進(jìn)一步的優(yōu)化是向量化加載。CUDA支持float4、int4等向量類型一次加載16字節(jié)。這樣每個(gè)線程處理4個(gè)floatwarp一次訪問(wèn)512字節(jié)transaction數(shù)量不變但指令數(shù)減少適合memory-bound的kernel。__global__ void vectorizedCopy(float4* out, float4* in, int n4) { int idx threadIdx.x blockIdx.x * blockDim.x; if (idx n4) { out[idx] in[idx]; } }注意使用float4要求數(shù)據(jù)地址16字節(jié)對(duì)齊。如果你用cudaMalloc分配內(nèi)存默認(rèn)就是256字節(jié)對(duì)齊沒(méi)問(wèn)題。但如果是偏移過(guò)的指針就要小心。實(shí)測(cè)數(shù)據(jù)在A100上標(biāo)量copy的帶寬大約是1.2 TB/sfloat4 copy能到1.5 TB/s左右提升了25%。這個(gè)差距在大規(guī)模數(shù)據(jù)處理里非?？捎^。3.2 Shared Memory Tiling矩陣乘法的經(jīng)典優(yōu)化矩陣乘法是CUDA優(yōu)化的“Hello World”。樸素版本的矩陣乘法每個(gè)線程計(jì)算C的一個(gè)元素需要讀取A的一行和B的一列。假設(shè)矩陣是N×N那總讀取量是2N3而計(jì)算量也是N3算術(shù)強(qiáng)度只有0.5完全是memory-bound。Tiling的思路是把A和B分成小塊每個(gè)block負(fù)責(zé)計(jì)算C的一個(gè)tile。block先把A和B對(duì)應(yīng)的tile加載到shared memory然后每個(gè)線程從shared memory里讀取數(shù)據(jù)做乘加。這樣global memory的讀取量降到2N3/tile_size算術(shù)強(qiáng)度提升了tile_size倍。#define TILE 16 __global__ void matmulTiled(float* C, float* A, float* B, int N) { __shared__ float As[TILE][TILE]; __shared__ float Bs[TILE][TILE]; int bx blockIdx.x, by blockIdx.y; int tx threadIdx.x, ty threadIdx.y; float sum 0.0f; for (int t 0; t N / TILE; t) { As[ty][tx] A[(by * TILE ty) * N t * TILE tx]; Bs[ty][tx] B[(t * TILE ty) * N bx * TILE tx]; __syncthreads(); for (int k 0; k TILE; k) { sum As[ty][k] * Bs[k][tx]; } __syncthreads(); } C[(by * TILE ty) * N bx * TILE tx] sum; }這段代碼里有兩個(gè)關(guān)鍵點(diǎn)一是加載到shared memory時(shí)global memory訪問(wèn)是合并的二是__syncthreads()保證所有線程都加載完了再計(jì)算計(jì)算完了再加載下一輪。實(shí)測(cè)1024×1024的矩陣乘法樸素版本大約2.5mstiled版本TILE16降到0.8msTILE32能到0.5ms左右。再往上調(diào)TILEshared memory用量增加occupancy下降收益遞減。3.3 Bank Conflictshared memory的隱形陷阱Shared memory被分成32個(gè)bank每個(gè)bank寬度4字節(jié)。如果warp里的32個(gè)線程訪問(wèn)的地址落在不同的bank上那一次就能完成如果多個(gè)線程訪問(wèn)同一個(gè)bank的不同地址就會(huì)發(fā)生bank conflict需要串行處理。最典型的情況是訪問(wèn)As[ty][k]這種按列訪問(wèn)。假設(shè)TILE32As是32×32的float數(shù)組那么As[0][0]和As[1][0]的地址相差32個(gè)float即128字節(jié)正好是32個(gè)bank的整數(shù)倍所以它們落在同一個(gè)bank上。如果warp里的線程同時(shí)訪問(wèn)As[0][0]、As[1][0]、...、As[31][0]那就是32路bank conflict性能直接崩掉。解決辦法是padding把As聲明為_(kāi)_shared__ float As[TILE][TILE1]多出一列。這樣As[0][0]和As[1][0]的地址相差33個(gè)float不再對(duì)齊到bank邊界conflict就消失了。__shared__ float As[TILE][TILE 1]; __shared__ float Bs[TILE][TILE 1];這個(gè)技巧簡(jiǎn)單但極其有效。我在一個(gè)圖像卷積的kernel里加了padding之后shared memory的吞吐量提升了將近3倍。3.4 減少分支分歧讓warp走同一條路前面提到warp divergence。如果你的kernel里有大量if-else而且條件依賴于threadIdx那warp里的線程就會(huì)走不同路徑串行執(zhí)行。解決辦法有幾種一是把條件改成warp-level的判斷讓整個(gè)warp走同一條路。比如if (threadIdx.x / 32 0)而不是if (threadIdx.x % 2 0)。二是用predication謂詞執(zhí)行代替分支。簡(jiǎn)單的if-else如果兩個(gè)分支都很短編譯器可能會(huì)自動(dòng)轉(zhuǎn)成predication但復(fù)雜分支就不行了。三是重新組織數(shù)據(jù)布局讓同一個(gè)warp的線程處理相同類型的數(shù)據(jù)。比如在稀疏矩陣計(jì)算里把非零元素按行分組每個(gè)warp處理一行就能減少分歧。實(shí)測(cè)在一個(gè)帶條件判斷的粒子模擬kernel里通過(guò)重新排序粒子讓同一warp的粒子處于相同狀態(tài)性能提升了40%。4. 實(shí)操?gòu)牧銉?yōu)化一個(gè)矩陣轉(zhuǎn)置kernel4.1 樸素版本為什么它慢得離譜矩陣轉(zhuǎn)置看起來(lái)很簡(jiǎn)單out[j][i] in[i][j]。但就是這個(gè)簡(jiǎn)單的操作樸素實(shí)現(xiàn)能慢到讓你懷疑顯卡壞了。__global__ void transposeNaive(float* out, float* in, int N) { int i blockIdx.y * blockDim.y threadIdx.y; int j blockIdx.x * blockDim.x threadIdx.x; if (i N j N) { out[j * N i] in[i * N j]; } }問(wèn)題在于讀取in的時(shí)候warp里的線程訪問(wèn)的是連續(xù)地址j連續(xù)這是合并的。但寫入out的時(shí)候地址是j*Niwarp里的線程j不同地址間隔是N完全無(wú)法合并。每個(gè)線程的寫入都會(huì)產(chǎn)生一個(gè)獨(dú)立的transaction帶寬利用率極低。在1024×1024的矩陣上這個(gè)kernel的耗時(shí)大約是1.8ms有效帶寬只有不到200 GB/s而A100的峰值帶寬是1.5 TB/s以上。4.2 Tiled轉(zhuǎn)置用shared memory做中轉(zhuǎn)優(yōu)化的思路是先用合并的方式把數(shù)據(jù)讀進(jìn)shared memory再?gòu)膕hared memory里按轉(zhuǎn)置后的順序讀出來(lái)合并地寫入global memory。#define TILE 32 __global__ void transposeTiled(float* out, float* in, int N) { __shared__ float tile[TILE][TILE 1]; int x blockIdx.x * TILE threadIdx.x; int y blockIdx.y * TILE threadIdx.y; if (x N y N) { tile[threadIdx.y][threadIdx.x] in[y * N x]; } __syncthreads(); x blockIdx.y * TILE threadIdx.x; y blockIdx.x * TILE threadIdx.y; if (x N y N) { out[y * N x] tile[threadIdx.x][threadIdx.y]; } }注意幾個(gè)細(xì)節(jié)一是shared memory數(shù)組加了paddingTILE1避免讀取tile[threadIdx.x][threadIdx.y]時(shí)的bank conflict二是__syncthreads()保證所有數(shù)據(jù)都加載完了再寫出三是寫出時(shí)的索引交換了blockIdx.x和blockIdx.y保證寫入也是合并的。實(shí)測(cè)同樣的1024×1024矩陣tiled版本的耗時(shí)降到0.3ms左右有效帶寬超過(guò)1 TB/s提升了將近6倍。4.3 參數(shù)調(diào)優(yōu)block size和grid size怎么選TILE32意味著block size是32×321024個(gè)線程。這是CUDA允許的最大block size。用滿1024的好處是shared memory利用率高但缺點(diǎn)是occupancy可能受限——每個(gè)SM最多支持2048個(gè)線程1024的block只能放2個(gè)如果寄存器用量高可能只能放1個(gè)。我試過(guò)TILE16block size 256和TILE32block size 1024在A100上TILE32略快但在一些較老的顯卡上TILE16反而更好因?yàn)閛ccupancy更高。所以沒(méi)有絕對(duì)的最優(yōu)值需要根據(jù)目標(biāo)硬件實(shí)測(cè)。Grid size的計(jì)算對(duì)于N×N矩陣gridDim.x (N TILE - 1) / TILEgridDim.y同理。注意邊界處理當(dāng)N不是TILE的整數(shù)倍時(shí)需要加if判斷。另外一個(gè)小技巧如果矩陣很大可以考慮用cudaMemcpy2D或者cudaMallocPitch來(lái)分配帶padding的行這樣即使不用shared memory也能在一定程度上改善訪問(wèn)模式。但shared memory tiling仍然是更通用的方案。5. 性能分析與調(diào)試用數(shù)據(jù)說(shuō)話5.1 Nsight Compute定位瓶頸的利器寫完kernel之后不要憑感覺(jué)猜哪里慢。用Nsight Compute跑一下它會(huì)告訴你kernel是compute-bound還是memory-boundoccupancy是多少有沒(méi)有bank conflictwarp divergence嚴(yán)重不嚴(yán)重。常用的幾個(gè)指標(biāo)sm__throughput.avg.pct_of_peak_sustained_elapsedSM整體吞吐量gpu__compute_memory_throughput.avg.pct_of_peak_sustained_elapsed內(nèi)存吞吐量l1tex__data_bank_conflicts_pipe_lsu_mem_shared_op_ld.sumshared memory bank conflict次數(shù)smsp__thread_inst_executed_per_inst_executed.ratio平均每個(gè)指令執(zhí)行的線程數(shù)低于32說(shuō)明有divergence如果內(nèi)存吞吐量接近100%而SM吞吐量很低那就是memory-bound優(yōu)化方向是減少內(nèi)存訪問(wèn)、提高合并度。反過(guò)來(lái)就是compute-bound優(yōu)化方向是減少指令數(shù)、用更快的指令比如用fma代替乘加分離。5.2 常見(jiàn)性能陷阱速查表癥狀可能原因排查方法解決思路帶寬利用率低于50%訪問(wèn)不合并檢查warp內(nèi)地址是否連續(xù)調(diào)整索引順序或用shared memoryshared memory吞吐量低bank conflict看bank conflict計(jì)數(shù)加padding或調(diào)整訪問(wèn)模式occupancy低但性能不差寄存器用量高看register per thread不一定要降看是否memory-bound性能隨數(shù)據(jù)量增大急劇下降緩存命中率低看L2命中率用tiling提高數(shù)據(jù)復(fù)用kernel啟動(dòng)開(kāi)銷大數(shù)據(jù)量太小看kernel執(zhí)行時(shí)間合并小kernel或改用CPU5.3 實(shí)測(cè)數(shù)據(jù)對(duì)比優(yōu)化前后的差距我在A100上跑了一組測(cè)試矩陣大小4096×4096float類型版本耗時(shí)(ms)有效帶寬(GB/s)相對(duì)提升樸素轉(zhuǎn)置28.51801.0xTiled轉(zhuǎn)置(TILE16)6.28304.6xTiled轉(zhuǎn)置(TILE32)4.810705.9xTiledfloat44.112507.0xfloat4版本是把每個(gè)線程處理4個(gè)元素進(jìn)一步減少了指令數(shù)。但要注意float4要求TILE能被4整除且邊界處理更復(fù)雜。6. 進(jìn)階話題從單卡到多卡從計(jì)算到通信6.1 多GPU編程數(shù)據(jù)怎么分通信怎么藏當(dāng)單卡放不下模型或者數(shù)據(jù)量太大時(shí)就需要多卡。CUDA提供了peer-to-peer access允許一塊卡直接訪問(wèn)另一塊卡的顯存。但P2P的帶寬遠(yuǎn)低于本地顯存所以核心思路是盡量讓計(jì)算在本地完成只交換必要的數(shù)據(jù)。常見(jiàn)的模式是每塊卡負(fù)責(zé)一部分?jǐn)?shù)據(jù)計(jì)算完成后把結(jié)果匯總到一塊卡上做后續(xù)處理。通信可以用cudaMemcpyPeerAsync來(lái)異步執(zhí)行和計(jì)算重疊。如果通信量不大還可以用NVLink如果硬件支持帶寬比PCIe高一個(gè)數(shù)量級(jí)。一個(gè)實(shí)用的技巧是把通信放在單獨(dú)的stream里和計(jì)算stream并行。這樣當(dāng)計(jì)算在跑的時(shí)候通信也在進(jìn)行總時(shí)間取決于兩者中較長(zhǎng)的那個(gè)而不是兩者之和。6.2 動(dòng)態(tài)并行與CUDA Graph減少啟動(dòng)開(kāi)銷動(dòng)態(tài)并行允許在kernel里啟動(dòng)另一個(gè)kernel適合遞歸或者自適應(yīng)細(xì)分的算法。但動(dòng)態(tài)并行的啟動(dòng)開(kāi)銷比host端啟動(dòng)還大所以只適合kernel執(zhí)行時(shí)間遠(yuǎn)大于啟動(dòng)開(kāi)銷的場(chǎng)景。CUDA Graph是另一種減少啟動(dòng)開(kāi)銷的方式。它把一系列kernel啟動(dòng)和內(nèi)存拷貝錄制成一個(gè)圖然后一次性提交。對(duì)于小kernel密集的場(chǎng)景能顯著降低CPU端的開(kāi)銷。我試過(guò)一個(gè)包含100個(gè)小kernel的流水線用CUDA Graph之后總時(shí)間從8ms降到5ms效果很明顯。6.3 混合精度與Tensor Core什么時(shí)候值得用如果你的計(jì)算涉及矩陣乘法而且能容忍一定的精度損失那Tensor Core是必須考慮的。A100的Tensor Core在FP16下的吞吐量是FP32的8倍以上。但Tensor Core對(duì)數(shù)據(jù)布局有要求需要把矩陣轉(zhuǎn)換成特定的格式比如行主序轉(zhuǎn)列主序或者用wmma API。混合精度的思路是用FP16做乘加用FP32做累加。這樣既利用了Tensor Core的高吞吐又保持了足夠的精度。實(shí)測(cè)在矩陣乘法上混合精度比純FP32快3-4倍精度損失在可接受范圍內(nèi)。但要注意不是所有算法都能容忍精度損失。比如涉及迭代求解的算法FP16的舍入誤差可能會(huì)累積放大。這時(shí)候要么用FP32要么用Kahan求和之類的補(bǔ)償技巧。7. 一些踩過(guò)的坑和實(shí)用建議第一個(gè)坑以為block size越大越好。實(shí)際上block size超過(guò)256之后收益往往遞減而且occupancy可能下降。我一般從128或256開(kāi)始試根據(jù)nsight的結(jié)果再調(diào)。第二個(gè)坑忘了__syncthreads()。shared memory的讀寫如果沒(méi)有同步會(huì)出現(xiàn)race condition結(jié)果時(shí)對(duì)時(shí)錯(cuò)非常難調(diào)試。記住寫shared memory之后、讀之前一定要同步。第三個(gè)坑在循環(huán)里頻繁調(diào)用cudaMalloc和cudaFree。這兩個(gè)操作非常慢應(yīng)該一次性分配好循環(huán)里復(fù)用。如果數(shù)據(jù)量動(dòng)態(tài)變化可以用內(nèi)存池或者cudaMallocAsync。第四個(gè)坑忽略錯(cuò)誤檢查。CUDA的API調(diào)用和kernel啟動(dòng)都可能出錯(cuò)但默認(rèn)不報(bào)錯(cuò)。每次調(diào)用后加cudaGetLastError()能省很多調(diào)試時(shí)間。第五個(gè)坑在debug模式下測(cè)性能。debug模式-G會(huì)禁用很多優(yōu)化性能數(shù)據(jù)沒(méi)有參考價(jià)值。測(cè)性能一定要用release模式-O3。最后分享一個(gè)小技巧如果你不確定某個(gè)優(yōu)化有沒(méi)有效果就寫兩個(gè)版本用cudaEvent計(jì)時(shí)跑100次取平均。不要憑感覺(jué)數(shù)據(jù)不會(huì)騙人。我見(jiàn)過(guò)太多人憑直覺(jué)優(yōu)化了半天結(jié)果性能反而下降了就是因?yàn)闆](méi)有做A/B測(cè)試。GPU編程的優(yōu)化空間很大但也沒(méi)有銀彈。理解硬件架構(gòu)、用工具定位瓶頸、小步快跑地驗(yàn)證這三條是我覺(jué)得最靠譜的路徑。希望這些經(jīng)驗(yàn)?zāi)軒湍闵僮咭恍澛贰?