CUDA系列學習(二)CUDA memory variables
本文來介紹CUDA的memory和變數存放,分為以下章節:
(一)、CPU Memory 結構
(二)、GPU Memory結構
(三)、CUDA Context
(四)、kernel設計
(五)、變數 & Memory
5.1 global arrays
5.2 global variables
5.3 Constant variables
5.4 Register
5.5 Local Array
5.6 Shared Memory
5.7 Texture Memory
5.8 總結
(一)、CPU Memory 結構
CPU提速主要依靠區域性性原理,即時間區域性性和空間區域性性。我們先看一下CPU的記憶體結構:
Data Access
先複習一下資料在這幾級儲存中的傳輸。作為資料transfer的基本單位,cache line的典型大小為8*8(8個變數,每個8bytes)=64bytes. 當一個cache想要load資料到暫存器時,檢查cache中的line,如果hit了就get到資料,否則將整條line從主存中去出來,(通常通過LRU)替換cache中一條line。暫存器傳資料到cache也一樣的過程。
Importance of Locality
上圖中可見在CPU中memory<--->L3 Cache傳輸頻寬為20GB/s, 除以64bytes/line得到傳輸記錄速度約300M line/s,約為300M*8= 2.4G double/s. 一般地,浮點數操作需要兩個輸入+1個輸出,那麼loading 3個數(3 lines)的代價為 100Mflops。如果一個line中的全部8個variables都被用到,那麼每秒浮點操作可以達到800Mflops。而CPU工作站典型為10 Gflops。這就要靠時間區域性性來重用資料了。
(二)、GPU Memory結構
Data Access
- Kepler GPU的cache line通常為128bytes(32個float or 16個double)。
- 資料傳輸頻寬最高250GB/s
- SMX的L2 cache統一1.5MB,L1 cache / shared memory有64KB
- 沒有CPU中的全域性快取一致性,所以幾乎沒有兩塊block更新相同的全域性陣列元素。
Importance of Locality
GPU對浮點數的操作速度可達1Tflops。和上面CPU的計算類似,GPU中memory<--->L2Cache傳輸頻寬為250GB/s, 除以128bytes/line得到傳輸記錄速度約2G line/s,約為2G*16= 32G double/s. 一般地,浮點數操作需要兩個輸入+1個輸出,那麼loading 3個數(3 lines)的代價為 670Mflops。如果一個line中的全部16個variables都被用到,那麼每秒浮點操作可以達到11Gflops。
這樣的話每進行一次資料到device的傳輸需要45flops(45次浮點操作)才能達到500Gflops. 所以很多演算法基本上不是卡在計算瓶頸,而是傳輸頻寬。
(三)、CUDA Context
一個CUDA Context類似於一個CPU程序。程式在Initialization的時候,runtime給每個device建立一個CUDA context,這個context在所有host threads中共享。driver API中的所有資源和action都封裝在一個CUDA context中,context被銷燬的時候系統自動清空這些資源,每個context擁有其自己的地址空間。所以,CUdeviceptr的value在不同context中會指向不同的記憶體空間。
一個host thread同一時刻只能用一個device context,每個host thread都有一個儲存當前contexts的stack。當一個context被cuCtxCreate()建立時,這個新的context被壓入棧(在棧頂),呼叫cuCtxPopCurrent() 可將這個context彈出來,然後這個context就會“漂”到其他host thread中再被壓入棧。
每個context都會維護一個count,表示有多少個threads在用。cuDtrCreate()令count = 1, cuCtxAttach()令count++,cuCtxDetach()令count--,cuCtxDestroy()令count = 0;一旦count=0,這個context就被銷燬。
(四)、kernel設計
我們在CUDA系列學習(一)中提到了GPU用的是SIMT cores,現在看一下它是如何進行執行緒管理的。每個SMX 多處理器在建立,管理,排程,執行的時候將threads每32個組成一組,稱為“wraps”。具體地,一個多處理器分配到多個blocks去執行的時候,它將blocks中的threads 分成wraps而且每個warp被一個warp scheduler來排程執行。一個warp一次執行一條相同指令,所以warp中所有threads同步執行是最有效的。那麼如果warp中的部分threads走上了資料相關的條件分支,warp就連續在各個branch上執行,暫停沒進入branch的threads。直到所有branch上的threads都執行完再合併了一起向下走。所以實現效能提升要注意儘量使warp內執行緒不要出現divergence。另外,注意這個branch divergence 之發生在warp內部;不同warp之間是獨立執行的。
看兩個kernel設計:
__global__ void kernel_1(float* x){ int tid = threadIdx.x + blockDim.x * blockIdx.x; x[tid] = threadIdx.x;}__global__ void kernel_2(float* x){ int tid = threadIdx.x + blockDim.x * blockIdx.x; x[1000*tid] = threadIdx.x;}
kernel_1中一個warp的32個thread訪問x的相鄰元素,即x[0]~x[31]在相同的cache line, 就是一個好的transfer;
kernel_2中訪問不連續記憶體,就要請求不同cache line,嚴重影響performance
(五)、變數 & Memory
上一篇中我們提到了memory由host memory和device memory組成,每部分尤其自己獨立的記憶體空間。Kernel跑在device memory上,所以runtime提供了分配,釋放,複製 device memory 和device <-->host 間transfer data的函式。
5.1 global arrays
global arrays:
- 儲存在/佔用device memory
- 由host code(非kernel部分code)宣告
- 一直存在,直到被host code釋放
- 因為所有block執行順序不定,所以如果一個block修改了一個數組元素,其他block就不能再對該元素進行讀寫
5.2 global variables
宣告前加識別符號__device__,表示變數要放在device上了 e.g. __device__ int reduction_lock=0;
__shared__(見4.6)和__constant__(見4.3)中至多有一個跟在__device__後面同時使用,標明用哪塊memory空間,如果這兩個都沒寫,則:
- 變數可以被grid內的所有threads讀寫
- 與application同生死
- 也可以定義為array,但是必須指定size
- 可以在host code中通過以下函式讀寫:
1. cudaMemcpyToSymbol;
2. cudaMemcpyFromSymbol;
3. cudaMemcpy + cudaGetSymbolAddress
Demo Code:
// float scalar__device__ float devData;float value = 3.14f;cudaMemcpyToSymbol(devData, &value, sizeof(float));//cudaMemcpyToSymbol(const char* symbol, const void* src, size_t count, size_t offset = 0, enum cudaMemcpyKind)// float array__device__ float* devPointer;float* ptr;cudaMalloc(&ptr, 256 * sizeof(float));cudaMemcpyToSymbol(devPointer, &ptr, sizeof(ptr));
5.3 Constant variables <常用>
- 哪裡宣告隨便,宣告前加識別符號__constant__
- 與application同生死
- grid內所有thread可直接讀(不可update),在host code中通過以下函式初始化
1. cudaMemcpyToSymbol;
2. cudaMemcpyFromSymbol;
3. cudaMemcpy + cudaGetSymbolAddress
Demo Code:
__constant__ float constData[256];float data[256];cudaMemcpyToSymbol(constData, data, sizeof(data)); //cudaMemcpyToSymbol(const char* symbol, const void* src, size_t count, size_t offset = 0, enum cudaMemcpyKind)cudaMemcpyFromSymbol(data, constData, sizeof(data)); //cudaMemcpyFromSymbol(const char* dst, const void* src_symbol, size_t count, size_t offset = 0, enum cudaMemcpyKind)
5.4 Register
- 預設一個kernel中的所有內部變數都存在register中
- 64K 32-bit registers per SMX
- up to 63 registers per thread (up to 255 for K20 / K40)
這時有64K/63 = 1024個threads (256個threads for K20 / K40)
- up to 2048 threads (at most 1024 per thread block)
這時每個thread有32個register
- not much difference between “fat” and “thin” threads
- 如果程式需要更多的register呢?就“spill over”到L1 cache,這樣訪問速度就慢了,我們要儘量避免spill
5.5 Local Array
指kernel code中宣告的陣列。
- 簡單情況下,編譯器會將小陣列float a[3]轉換成3個標量registers:a0,a1,a2作處理
- 複雜的情況,會將array放到L1(16KB),只能放4096個32-bit的變數,如果有1024個執行緒,每個執行緒只能分配放4個變數。
5.6 Shared Memory
前面加識別符號__shared__ e.g. __shared__ int x_dim;
- 要佔用thread block的shared memory space.
- 要比global memory快很多,所以只要有機會就把global memory整成shared memory
- 與block同生死
- thread block內所有threads共用(可讀可寫)
- 啥時侯用呢?當所有threads訪問都是同一個值的時候,這樣就避免用register了
但是有問題就是,如果一個thread block有多個warp(上一篇blog中提到的概念,block中的thread每32個被分到一個warp,最後一個不足32個thread也沒關係,同樣形成一個warp),各warp執行指令順序是不定的,那麼久需要執行緒同步機制,用指令__syncthreads(); 插入一個“barrier”,所有wrap執行到這個barrier之前沒有thread/warp能夠越過去。
Kepler GPU給L1 Cache + shared memory總共64KB,可以分為16+48,32+32,48+16;這個split可以通過cudaFuncSetCacheConfig()或cudaDeviceSetCacheConfig()設定,預設給shared memroy 48KB。這個具體情況看程式了。
下面通過一個經典例子來看shared memory作用:矩陣乘法
目的:實現C=A*B,方法:c[i,j] = A[i,:] * B[:,j],
其中矩陣用row-major表示,即c[i,j] = *(c.elements + i*c.width + j)
1. 不用shared memory優化版:
設A為m*t的矩陣;B為t*n的矩陣;
每個執行緒讀取A的一行,B的一列,計算C的對應值;
所以這樣需要從global memory中讀n次A,m次B。
// Matrices are stored in row-major order:// M(row, col) = *(M.elements + row * M.width + col)typedef struct { int width; int height; float* elements;} Matrix;// Thread block size#define BLOCK_SIZE 16// Forward declaration of the matrix multiplication kernel__global__ void MatMulKernel(const Matrix, const Matrix, Matrix);// Matrix multiplication - Host code// Matrix dimensions are assumed to be multiples of BLOCK_SIZEvoid MatMul(const Matrix A, const Matrix B, Matrix C){ // Load A and B to device memory Matrix d_A; d_A.width = A.width; d_A.height = A.height; size_t size = A.width * A.height * sizeof(float); cudaMalloc(&d_A.elements, size); cudaMemcpy(d_A.elements, A.elements, size, cudaMemcpyHostToDevice); Matrix d_B; d_B.width = B.width; d_B.height = B.height; size = B.width * B.height * sizeof(float); cudaMalloc(&d_B.elements, size); cudaMemcpy(d_B.elements, B.elements, size, cudaMemcpyHostToDevice); // Allocate C in device memory Matrix d_C; d_C.width = C.width; d_C.height = C.height; size = C.width * C.height * sizeof(float); cudaMalloc(&d_C.elements, size); // Invoke kernel dim3 dimBlock(BLOCK_SIZE, BLOCK_SIZE); dim3 dimGrid(B.width / dimBlock.x, A.height / dimBlock.y); MatMulKernel<<<dimGrid, dimBlock>>>(d_A, d_B, d_C); // Read C from device memory cudaMemcpy(C.elements, Cd.elements, size, cudaMemcpyDeviceToHost); } // Free device memory cudaFree(d_A.elements); cudaFree(d_B.elements); cudaFree(d_C.elements);}// Matrix multiplication kernel called by MatMul()__global__ void MatMulKernel(Matrix A, Matrix B, Matrix C){ // Each thread computes one element of C // by accumulating results into Cvalue float Cvalue = 0; int row = blockIdx.y * blockDim.y + threadIdx.y; int col = blockIdx.x * blockDim.x + threadIdx.x; for (int e = 0; e < A.width; ++e) Cvalue += A.elements[row * A.width + e]* B.elements[e * B.width + col]; C.elements[row * C.width + col] = Cvalue;}
2. 利用shared memory
每個thread block負責計算一個子矩陣Csub, 其中每個thread負責計算Csub中的一個元素。如下圖所示。為了將fit裝置資源,A,B都分割成很多block_size維的方形matrix,Csub將這些方形matrix的乘積求和而得。每次計算一個乘積時,先將兩個對應方形矩陣從global memory 載入 shared memory(一個thread負責載入A, B兩個sub matrix的元素),然後每個thread計算乘積的一個元素,再由每個thread將這些product加和,存入一個register,最後一次性寫入global memory。計算時注意同步,詳見程式碼。
設A為m*t的矩陣;B為t*n的矩陣;
這樣呢,A只從global memory讀了n/block_size次,B只讀了m/block_size次;
Kernel Code:
__global__ void MatMulKernel(Matrix A, Matrix B, Matrix C){ // Block row and column int blockRow = blockIdx.y; int blockCol = blockIdx.x; // Each thread block computes one sub-matrix Csub of C Matrix Csub = GetSubMatrix(C, blockRow, blockCol); // Each thread computes one element of Csub by accumulating results into Cvalue float Cvalue = 0; // Thread row and column within Csub int row = threadIdx.y; int col = threadIdx.x; // Loop over all the sub-matrices of A and B that are // required to compute Csub // Multiply each pair of sub-matrices together // and accumulate the results for (int m = 0; m < (A.width / BLOCK_SIZE); ++m) { // Get sub-matrix Asub of A Matrix Asub = GetSubMatrix(A, blockRow, m); // Get sub-matrix Bsub of B Matrix Bsub = GetSubMatrix(B, m, blockCol); // Shared memory used to store Asub and Bsub respectively __shared__ float As[BLOCK_SIZE][BLOCK_SIZE]; __shared__ float Bs[BLOCK_SIZE][BLOCK_SIZE]; // Load Asub and Bsub from device memory to shared memory // Each thread loads one element of each sub-matrix As[row][col] = GetElement(Asub, row, col); Bs[row][col] = GetElement(Bsub, row, col); // Synchronize to make sure the sub-matrices are loaded // before starting the computation __syncthreads(); // Multiply Asub and Bsub together for (int e = 0; e < BLOCK_SIZE; ++e) Cvalue += As[row][e] * Bs[e][col]; // Synchronize to make sure that the preceding // computation is done before loading two new // sub-matrices of A and B in the next iteration __syncthreads(); } // Write Csub to device memory // Each thread writes one element SetElement(Csub, row, col, Cvalue);}
Host Code:
// Invoke kerneldim3 dimBlock(BLOCK_SIZE, BLOCK_SIZE);dim3 dimGrid(B.width / dimBlock.x, A.height / dimBlock.y);MatMulKernel<<<dimGrid, dimBlock>>>(d_A, d_B, d_C);
5.7 Texture memory
前面加識別符號const __restrict__, 之所以叫texture是因為之前用texture memory想服務於純graphics的應用。
不同於shared memory,對texture memory, 不同執行緒可以訪問到不同value。K20/K40中texture cache有48KB。
5.8 總結
綜上,每個block內有以下資源:
- threads
- registers (registers per thread * number of threads)
- shared memory
這些決定了一個SMX上能同時執行多少個blocks(最多16個)。
參考:
6. Fermi 架構白皮書(GPU繼承了Fermi的很多架構特點)