cuda入门——改良第一个 CUDA程序
在上篇中，我們做了一個計算一大堆數字的平方和的程序。不過，我們也提到這個程序的執行效率并不理想。當然，實際上來說，如果只是要做計算平方和的動作，用 CPU 做會比用 GPU 快得多。這是因為平方和的計算并不需要太多運算能力，所以幾乎都是被內存帶寬所限制。因此，光是把數據復制到顯卡內存上的這個動作，所需要的時間，可能已經和直接在 CPU 上進行計算差不多了。
不過，如果進行平方和的計算，只是一個更復雜的計算過程的一部份的話，那么當然在 GPU 上計算還是有它的好處的。而且，如果數據已經在顯卡內存上（例如在 GPU 上透過某種算法產生），那么，使用 GPU 進行這樣的運算，還是會比較快的。
剛才也提到了，由于這個計算的主要瓶頸是內存帶寬，所以，理論上顯卡的內存帶寬是相當大的。這里我們就來看看，倒底我們的第一個程序，能利用到多少內存帶寬。
| 程序的并行化 |
我們的第一個程序，并沒有利用到任何并行化的功能。整個程序只有一個 thread。在 GeForce 8800GT 上面，在 GPU 上執行的部份（稱為 "kernel"）大約花費 640M 個頻率。GeForce 8800GT 的執行單元的頻率是 1.5GHz，因此這表示它花費了約 0.43 秒的時間。1M 個 32 bits 數字的數據量是 4MB，因此，這個程序實際上使用的內存帶寬，只有 9.3MB/s 左右！這是非常糟糕的表現。
為什么會有這樣差的表現呢？這是因為 GPU 的架構特性所造成的。在 CUDA 中，一般的數據復制到的顯卡內存的部份，稱為 global memory。這些內存是沒有 cache 的，而且，存取 global memory 所需要的時間（即 latency）是非常長的，通常是數百個 cycles。由于我們的程序只有一個 thread，所以每次它讀取 global memory 的內容，就要等到實際讀取到數據、累加到 sum 之后，才能進行下一步。這就是為什么它的表現會這么的差。
由于 global memory 并沒有 cache，所以要避開巨大的 latency 的方法，就是要利用大量的 threads。假設現在有大量的 threads 在同時執行，那么當一個 thread 讀取內存，開始等待結果的時候，GPU 就可以立刻切換到下一個 thread，并讀取下一個內存位置。因此，理想上當 thread 的數目夠多的時候，就可以完全把 global memory 的巨大 latency 隱藏起來了。
要怎么把計算平方和的程序并行化呢？最簡單的方法，似乎就是把數字分成若干組，把各組數字分別計算平方和后，最后再把每組的和加總起來就可以了。一開始，我們可以把最后加總的動作，由 CPU 來進行。
首先，在 first_cuda.cu 中，在 #define DATA_SIZE 的后面增加一個 #define，設定 thread 的數目：
#define DATA_SIZE??? 1048576
#define THREAD_NUM?? 256
接著，把 kernel 程序改成：
__global__ static void sumOfSquares(int *num, int* result,
??? clock_t* time)
{
??? const int tid = threadIdx.x;
??? const int size = DATA_SIZE / THREAD_NUM;
??? int sum = 0;
??? int i;
??? clock_t start;
??? if(tid == 0) start = clock();
??? for(i = tid * size; i < (tid + 1) * size; i++) {
?????? sum += num[i] * num[i];
??? }
??? result[tid] = sum;
??? if(tid == 0) *time = clock() - start;
}
程序里的 threadIdx 是 CUDA 的一個內建的變量，表示目前的 thread 是第幾個 thread（由 0 開始計算）。以我們的例子來說，會有 256 個 threads，所以同時會有 256 個 sumOfSquares 函式在執行，但每一個的 threadIdx.x 則分別會是 0 ~ 255。利用這個變量，我們就可以讓每一份函式執行時，對整個數據不同的部份計算平方和。另外，我們也讓計算時間的動作，只在 thread 0（即 threadIdx.x = 0 的時候）進行。
同樣的，由于會有 256 個計算結果，所以原來存放 result 的內存位置也要擴大。把 main 函式中的中間部份改成：
??? int* gpudata, *result;
??? clock_t* time;
??? cudaMalloc((void**) &gpudata, sizeof(int) * DATA_SIZE);
??? cudaMalloc((void**) &result, sizeof(int) * THREAD_NUM);
??? cudaMalloc((void**) &time, sizeof(clock_t));
??? cudaMemcpy(gpudata, data, sizeof(int) * DATA_SIZE,
??? ??? cudaMemcpyHostToDevice);
??? sumOfSquares<<<1, THREAD_NUM, 0>>>(gpudata, result, time);
??? int sum[THREAD_NUM];
??? clock_t time_used;
??? cudaMemcpy(&sum, result, sizeof(int) * THREAD_NUM,
??? ??? cudaMemcpyDeviceToHost);
??? cudaMemcpy(&time_used, time, sizeof(clock_t),
??? ??? cudaMemcpyDeviceToHost);
??? cudaFree(gpudata);
??? cudaFree(result);
??? cudaFree(time);
可以注意到我們在呼叫 sumOfSquares 函式時，指定 THREAD_NUM 為 thread 的數目。最后，在 CPU 端把計算好的各組數據的平方和進行加總：
??? int final_sum = 0;
??? for(int i = 0; i < THREAD_NUM; i++) {
??? ??? final_sum += sum[i];
??? }
??? printf("sum: %d time: %d\n", final_sum, time_used);
??? final_sum = 0;
??? for(int i = 0; i < DATA_SIZE; i++) {
??????? sum += data[i] * data[i];
??? }
??? printf("sum (CPU): %d\n", final_sum);
編譯后執行，確認結果和原來相同。
這個版本的程序，在 GeForce 8800GT 上執行，只需要約 8.3M cycles，比前一版程序快了 77 倍！這就是透過大量 thread 來隱藏 latency 所帶來的效果。
不過，如果計算一下它使用的內存帶寬，就會發現其實仍不是很理想，大約只有 723MB/s 而已。這和 GeForce 8800GT 所具有的內存帶寬是很大的差距。為什么會這樣呢？
| 內存的存取模式 |
顯卡上的內存是 DRAM，因此最有效率的存取方式，是以連續的方式存取。前面的程序，雖然看起來是連續存取內存位置（每個 thread 對一塊連續的數字計算平方和），但是我們要考慮到實際上 thread 的執行方式。前面提過，當一個 thread 在等待內存的數據時，GPU 會切換到下一個 thread。也就是說，實際上執行的順序是類似
??? thread 0 -> thread 1 -> thread 2 -> ...
因此，在同一個 thread 中連續存取內存，在實際執行時反而不是連續了。要讓實際執行結果是連續的存取，我們應該要讓 thread 0 讀取第一個數字，thread 1 讀取第二個數字…依此類推。所以，我們可以把 kernel 程序改成如下：
__global__ static void sumOfSquares(int *num, int* result,
??? clock_t* time)
{
??? const int tid = threadIdx.x;
??? int sum = 0;
??? int i;
??? clock_t start;
??? if(tid == 0) start = clock();
??? for(i = tid; i < DATA_SIZE; i += THREAD_NUM) {
?????? sum += num[i] * num[i];
??? }
??? result[tid] = sum;
??? if(tid == 0) *time = clock() - start;
}
編譯后執行，確認結果相同。
僅僅是這樣簡單的修改，實際執行的效率就有很大的差別。在 GeForce 8800GT 上，上面的程序執行需要的頻率是 2.6M cycles，又比前一版程序快了三倍。不過，這樣仍只有 2.3GB/s 的帶寬而已。
這是因為我們使用的 thread 數目還是不夠多的原因。理論上 256 個 threads 最多只能隱藏 256 cycles 的 latency。但是 GPU 存取 global memory 時的 latency 可能高達 500 cycles 以上。如果增加 thread 數目，就可以看到更好的效率。例如，可以把 THREAD_NUM 改成 512。在 GeForce 8800GT 上，這可以讓執行花費的時間減少到 1.95M cycles。有些改進，但是仍不夠大。不幸的是，目前 GeForce 8800GT 一個 block 最多只能有 512 個 threads，所以不能再增加了，而且，如果 thread 數目增加太多，那么在 CPU 端要做的最后加總工作也會變多。
| 更多的并行化 |
前面提到了 block。在之前介紹呼叫 CUDA 函式時，也有提到 "block 數目" 這個參數。到目前為止，我們都只使用一個 block。究竟 block 是什么呢？
在 CUDA 中，thread 是可以分組的，也就是 block。一個 block 中的 thread，具有一個共享的 shared memory，也可以進行同步工作。不同 block 之間的 thread 則不行。在我們的程序中，其實不太需要進行 thread 的同步動作，因此我們可以使用多個 block 來進一步增加 thread 的數目。
首先，在 #define DATA_SIZE 的地方，改成如下：
#define DATA_SIZE?? 1048576
#define BLOCK_NUM?? 32
#define THREAD_NUM?? 256
這表示我們會建立 32 個 blocks，每個 blocks 有 256 個 threads，總共有 32*256 = 8192 個 threads。
接著，我們把 kernel 部份改成：
__global__ static void sumOfSquares(int *num, int* result,
??? clock_t* time)
{
??? const int tid = threadIdx.x;
??? const int bid = blockIdx.x;
??? int sum = 0;
??? int i;
??? if(tid == 0) time[bid] = clock();
??? for(i = bid * THREAD_NUM + tid; i < DATA_SIZE;
??? ??? i += BLOCK_NUM * THREAD_NUM) {
?????? sum += num[i] * num[i];
??? }
??? result[bid * THREAD_NUM + tid] = sum;
??? if(tid == 0) time[bid + BLOCK_NUM] = clock();
}
blockIdx.x 和 threadIdx.x 一樣是 CUDA 內建的變量，它表示的是目前的 block 編號。另外，注意到我們把計算時間的方式改成每個 block 都會記錄開始時間及結束時間。
main 函式部份，修改成：
??? int* gpudata, *result;
??? clock_t* time;
??? cudaMalloc((void**) &gpudata, sizeof(int) * DATA_SIZE);
??? cudaMalloc((void**) &result,
??? ??? sizeof(int) * THREAD_NUM * BLOCK_NUM);
??? cudaMalloc((void**) &time, sizeof(clock_t) * BLOCK_NUM * 2);
??? cudaMemcpy(gpudata, data, sizeof(int) * DATA_SIZE,
??? ??? cudaMemcpyHostToDevice);
??? sumOfSquares<<<BLOCK_NUM, THREAD_NUM, 0>>>(gpudata, result,
??? ??? time);
??? int sum[THREAD_NUM * BLOCK_NUM];
??? clock_t time_used[BLOCK_NUM * 2];
??? cudaMemcpy(&sum, result, sizeof(int) * THREAD_NUM * BLOCK_NUM,
??? ??? cudaMemcpyDeviceToHost);
??? cudaMemcpy(&time_used, time, sizeof(clock_t) * BLOCK_NUM * 2,
??? ??? cudaMemcpyDeviceToHost);
??? cudaFree(gpudata);
??? cudaFree(result);
??? cudaFree(time);
??? int final_sum = 0;
??? for(int i = 0; i < THREAD_NUM * BLOCK_NUM; i++) {
??? ??? final_sum += sum[i];
??? }
??? clock_t min_start, max_end;
??? min_start = time_used[0];
??? max_end = time_used[BLOCK_NUM];
??? for(int i = 1; i < BLOCK_NUM; i++) {
??? ??? if(min_start > time_used[i])
??? ??? ??? min_start = time_used[i];
??? ??? if(max_end < time_used[i + BLOCK_NUM])
??? ??? ??? max_end = time_used[i + BLOCK_NUM];
??? }
??? printf("sum: %d time: %d\n", final_sum, max_end - min_start);
基本上我們只是把 result 的大小變大，并修改計算時間的方式，把每個 block 最早的開始時間，和最晚的結束時間相減，取得總運行時間。
這個版本的程序，執行的時間減少很多，在 GeForce 8800GT 上只需要約 150K cycles，相當于 40GB/s 左右的帶寬。不過，它在 CPU 上執行的部份，需要的時間加長了（因為 CPU 現在需要加總 8192 個數字）。為了避免這個問題，我們可以讓每個 block 把自己的每個 thread 的計算結果進行加總。
| Thread 的同步 |
前面提過，一個 block 內的 thread 可以有共享的內存，也可以進行同步。我們可以利用這一點，讓每個 block 內的所有 thread 把自己計算的結果加總起來。把 kernel 改成如下：
__global__ static void sumOfSquares(int *num, int* result,
??? clock_t* time)
{
??? extern __shared__ int shared[];
??? const int tid = threadIdx.x;
??? const int bid = blockIdx.x;
??? int i;
??? if(tid == 0) time[bid] = clock();
??? shared[tid] = 0;
??? for(i = bid * THREAD_NUM + tid; i < DATA_SIZE;
??? ??? i += BLOCK_NUM * THREAD_NUM) {
?????? shared[tid] += num[i] * num[i];
??? }
??? __syncthreads();
??? if(tid == 0) {
??? ??? for(i = 1; i < THREAD_NUM; i++) {
??? ??? ??? shared[0] += shared[i];
??? ??? }
??? ??? result[bid] = shared[0];
??? }
??? if(tid == 0) time[bid + BLOCK_NUM] = clock();
}
利用 __shared__ 聲明的變量表示這是 shared memory，是一個 block 中每個 thread 都共享的內存。它會使用在 GPU 上的內存，所以存取的速度相當快，不需要擔心 latency 的問題。
__syncthreads() 是一個 CUDA 的內部函數，表示 block 中所有的 thread 都要同步到這個點，才能繼續執行。在我們的例子中，由于之后要把所有 thread 計算的結果進行加總，所以我們需要確定每個 thread 都已經把結果寫到 shared[tid] 里面了。
接下來，把 main 函式的一部份改成：
??? int* gpudata, *result;
??? clock_t* time;
??? cudaMalloc((void**) &gpudata, sizeof(int) * DATA_SIZE);
??? cudaMalloc((void**) &result, sizeof(int) * BLOCK_NUM);
??? cudaMalloc((void**) &time, sizeof(clock_t) * BLOCK_NUM * 2);
??? cudaMemcpy(gpudata, data, sizeof(int) * DATA_SIZE,
??? ??? cudaMemcpyHostToDevice);
??? sumOfSquares<<<BLOCK_NUM, THREAD_NUM,
??? ??? THREAD_NUM * sizeof(int)>>>(gpudata, result, time);
??? int sum[BLOCK_NUM];
??? clock_t time_used[BLOCK_NUM * 2];
??? cudaMemcpy(&sum, result, sizeof(int) * BLOCK_NUM,
??? ??? cudaMemcpyDeviceToHost);
??? cudaMemcpy(&time_used, time, sizeof(clock_t) * BLOCK_NUM * 2,
??? ??? cudaMemcpyDeviceToHost);
??? cudaFree(gpudata);
??? cudaFree(result);
??? cudaFree(time);
??? int final_sum = 0;
??? for(int i = 0; i < BLOCK_NUM; i++) {
??? ??? final_sum += sum[i];
??? }
可以注意到，現在 CPU 只需要加總 BLOCK_NUM 也就是 32 個數字就可以了。
不過，這個程序由于在 GPU 上多做了一些動作，所以它的效率會比較差一些。在 GeForce 8800GT 上，它需要約 164K cycles。
當然，效率會變差的一個原因是，在這一版的程序中，最后加總的工作，只由每個 block 的 thread 0 來進行，但這并不是最有效率的方法。理論上，把 256 個數字加總的動作，是可以并行化的。最常見的方法，是透過樹狀的加法：
把 kernel 改成如下：
__global__ static void sumOfSquares(int *num, int* result,
??? clock_t* time)
{
??? extern __shared__ int shared[];
??? const int tid = threadIdx.x;
??? const int bid = blockIdx.x;
??? int i;
??? int offset = 1, mask = 1;
??? if(tid == 0) time[bid] = clock();
??? shared[tid] = 0;
??? for(i = bid * THREAD_NUM + tid; i < DATA_SIZE;
??? ??? i += BLOCK_NUM * THREAD_NUM) {
?????? shared[tid] += num[i] * num[i];
??? }
??? __syncthreads();
??? while(offset < THREAD_NUM) {
??? ??? if((tid & mask) == 0) {
??? ??? ??? shared[tid] += shared[tid + offset];
??? ??? }
??? ??? offset += offset;
??? ??? mask = offset + mask;
??? ??? __syncthreads();
??? }
??? if(tid == 0) {
??? ??? result[bid] = shared[0];
??? ??? time[bid + BLOCK_NUM] = clock();
??? }
}
后面的 while 循環就是進行樹狀加法。main 函式則不需要修改。
這一版的程序，在 GeForce 8800GT 上執行需要的時間，大約是 140K cycles（相當于約 43GB/s），比完全不在 GPU 上進行加總的版本還快！這是因為，在完全不在 GPU 上進行加總的版本，寫入到 global memory 的數據數量很大（8192 個數字），也對效率會有影響。所以，這一版程序不但在 CPU 上的運算需求降低，在 GPU 上也能跑的更快。
| 進一步改善 |
上一個版本的樹狀加法是一般的寫法，但是它在 GPU 上執行的時候，會有 share memory 的 bank conflict 的問題（詳情在后面介紹 GPU 架構時會提到）。采用下面的方法，可以避免這個問題：
??? offset = THREAD_NUM / 2;
??? while(offset > 0) {
??? ??? if(tid < offset) {
??? ??? ??? shared[tid] += shared[tid + offset];
??? ??? }
??? ??? offset >>= 1;
??? ??? __syncthreads();
??? }
這樣同時也省去了 mask 變數。因此，這個版本的執行的效率就可以再提高一些。在 GeForce 8800GT 上，這個版本執行的時間是約 137K cycles。當然，這時差別已經很小了。如果還要再提高效率，可以把樹狀加法整個展開：
??? if(tid < 128) { shared[tid] += shared[tid + 128]; }
??? __syncthreads();
??? if(tid < 64) { shared[tid] += shared[tid + 64]; }
??? __syncthreads();
??? if(tid < 32) { shared[tid] += shared[tid + 32]; }
??? __syncthreads();
??? if(tid < 16) { shared[tid] += shared[tid + 16]; }
??? __syncthreads();
??? if(tid < 8) { shared[tid] += shared[tid + 8]; }
??? __syncthreads();
??? if(tid < 4) { shared[tid] += shared[tid + 4]; }
??? __syncthreads();
??? if(tid < 2) { shared[tid] += shared[tid + 2]; }
??? __syncthreads();
??? if(tid < 1) { shared[tid] += shared[tid + 1]; }
??? __syncthreads();
當然這只適用于 THREAD_NUM 是 256 的情形。這樣可以再省下約 1000 cycles 左右（約 44GB/s）。
posted on 2009-04-28 19:53 GXW 閱讀(...) 評論(...) 編輯 收藏轉載于:https://www.cnblogs.com/Fancyboy2004/archive/2009/04/28/1445635.html
總結
以上是生活随笔為你收集整理的cuda入门——改良第一个 CUDA程序的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： 让ModalPopupExtender的
- 下一篇： 转载--让搜索引擎优化(SEO)更有效