[分享] CUDA 程式設計(7) -- 來玩泡泡吧

[a5000ml (咖啡裡的海洋藍)]

第七章  來玩泡泡吧

今天我們來玩泡泡排序法吧!

◆ 演算法 

        為了能夠平行化, 使用交錯的方式對配對的資料進行排序, 如下圖所示, 假設
        有 10 筆資料, 我們先對位址 (0 1),(2 3),... 的 0 based 資料對進行排序,
        然後再對位址 (1 2),(3 4),... 的 1 based 資料對進行排序, 如此一來,
        每次進程可讓最小的往前浮上兩個單位, 總共要執行 N/2 次才會浮到最前面,
        每次比對數量為 N-1, 所以是 O(N^2/2) 的泡泡法.

                (0  1)(2  3)(4  5)(6  7)(8  9)      0-based
                 0 (1  2)(3  4)(5  6)(7  8) 9       1-based

        若資料量為奇數時, 例如 N=11, 則執行的方式如下,

                (0  1)(2  3)(4  5)(6  7)(8  9) 10   0-based
                 0 (1  2)(3  4)(5  6)(7  8)(9  10)  1-based

        所以端點處和偶數時略有不同, 玩泡泡時要注意一下.


◆ 單一區塊範例

        以下示範用共享記憶體來做泡泡排序, 每個執行緒負責比較一對相鄰的資料,
        因為單一區塊最大的執行緒個數是 512, 所以這個版本最多的輸入資料數為

                N = 512*2+1 = 1025

        每次進程先比較 0-based 配對, 再比較 1-based 配對, 總共要執行 N/2 次,
        也就是 blockDim 次, 程式碼放在附錄 bubble1.cu, 在 GTX 8800 Ultra 上
        執行的結果如下 (host 為Intel Q6600)

                time[gpu]: 0.72 ms
                time[host]: 1.36 ms
                ratio: 1.88889 x
                ------------------------
                test[gpu]: pass
                test[host]: pass

        效能只有 1.8 倍, 因為它只用到一個區塊, 也就是只用到 G80 上面 16 個
        多處理器 (MP, multiprocessors) 的其中之一, 如果能夠使用很多個區塊
        不僅效能能提昇數倍, 而且也能加大排序資料量.


◆ 練習 

        (1) 將範例改成用全域記憶體的版本, 並和使用共享記憶體的範例比較效能.

        (2) 將泡泡排序改為多區塊版本, 首先將資料分段, 然後使用多區塊對每段
            進行排序, 其分段法仿照範例, 例如要排序 10*512 筆資料, 便可設定
            10 個區塊, 每個區塊使用 256 個執行緒做區塊的泡泡排序, 整個網格
            先做 0-based, 再做 256-based, 每次發出的區塊如下

                 0-based    (0~511),   (512~1023), (1024~1535), ...
                 256-based  (256~767), (768~1279), (1280~1791), ...

            因為每次最小的群組都會往前步進一個區塊, 所以對 10*512 而言, 只要
            讓它 loop 個 10 次, 就可以把資料全部排序好了.
            (Note: 網格發出可在 host 上控制, 結果於下禮拜公佈)

        (3) 試估計多區塊版本的效能, 並與結果比較. (需考慮 MP 個數)



========================================================================
bubble1.cu
========================================================================
#include <cuda.h>
#include <time.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

//----------------------------------------------
// 排序 N 個 float 元素 (N=1~1025)
// 使用到的記憶體大小為 SIZE 個 bytes
// 只使用單一區塊
// 區塊大小為 N/2
//----------------------------------------------
#define N      1024
#define SIZE   (N*sizeof(float))

#define GRID   1
#define BLOCK  (N/2)

#define testLoop 1000  //測試效能時的 loop 數

//----------------------------------------------
// 交換函式 (host 和 kernel 都可以使用) 
// 因為加了 __host__ 和 __device__ 兩個標籤
//----------------------------------------------
inline __host__ __device__ void swap(float& a, float& b){
        float c=a;
        a=b;
        b=c;
}

//----------------------------------------------
// 泡泡的 kernel (由小到大排列 N 個元素 a->r)
//----------------------------------------------
__global__ void bubble(float *r, float *a){
        //*** blockDim=N/2 ***
        int j=threadIdx.x;      //j=0,1,2,...blockDim-1 
        int k=2*threadIdx.x;    //k=0,2,4,...2*(blockDim-1) 配對的基底索引

        //配置共享記憶體
        __shared__ float s[N+20];

        //載入資料到共享記憶體
        __syncthreads();   //同步化執行緒, 加速載入速度 (合併讀取 coalesced)
        s[j]=a[j];         //使用全部執行緒一起載入前半段 (0~N/2-1)
        s[j+N/2]=a[j+N/2]; //使用全部執行緒一起載入後半段 (N/2~N-1)
        if(j==0){
                //若 N 為奇數時, 還要多載入一個尾巴, 只使用第 0 個執行緒
                s[N-1]=a[N-1];
        }

        //開始泡泡排序
        for(int loop=0; loop<=N/2; loop++){
                //排列 0 based 配對資料 (0,1) (2,3) (4,5) .... 
                __syncthreads();  //同步化確保共享記憶體已寫入
                if(s[k]>s[k+1]){
                        swap(s[k],s[k+1]);
                }

                //排列 1 based 配對資料 (1,2) (3,4) (5,6) .... 
                __syncthreads();  //同步化確保共享記憶體已寫入
                if(s[k+1]>s[k+2]){
                        if(k<N-2) //若 N 為偶數時, 最後一個執行緒不作用
                        swap(s[k+1],s[k+2]);
                }
        }

        //轉出資料到全域記憶體
        __syncthreads();
        r[j]=s[j];
        r[j+N/2]=s[j+N/2];
        if(j==0){
                r[N-1]=s[N-1];
        }

}


//----------------------------------------------
// 泡泡的 host 函數 
//----------------------------------------------
void bubble_host(float *r, float *a){
        //載入資料
        for(int k=0; k<N; k++){
                r[k]=a[k];
        }

        for(int loop=0; loop<=N/2; loop++){
                //排列 0 based 配對資料
                for(int k=0; k<N-1; k+=2){
                        if(r[k]>r[k+1]){
                                swap(r[k],r[k+1]);
                        }
                }
                //排列 1 based 配對資料
                for(int k=1; k<N-1; k+=2){
                        if(r[k]>r[k+1]){
                                swap(r[k],r[k+1]);
                        }
                }
        }
}


//----------------------------------------------
// 主程式
//----------------------------------------------
int main(){
        //配置 host 記憶體
        float *a=(float*)malloc(SIZE);
        float *b=(float*)malloc(SIZE);
        float *c=(float*)malloc(SIZE);


        //初始化
        for(int k=0; k<N; k++){
                a[k]=k;
                c[k]=0;
        }

        //對陣列 a 洗牌
        srand(time(0));
        for(int k=0; k<2*N; k++){
                int i=rand()%N;
                int j=rand()%N;
                swap(a[i],a[j]);
        }

        //配置 device 記憶體
        float  *ga, *gc;
        cudaMalloc((void**)&ga, SIZE);
        cudaMalloc((void**)&gc, SIZE);

        //載入 (順便載入 c 來清空裝置記憶體內容)
        cudaMemcpy(ga, a, SIZE, cudaMemcpyHostToDevice);
        cudaMemcpy(gc, c, SIZE, cudaMemcpyHostToDevice);


        //測試 kernel 效能
        double t0=(double)clock()/CLOCKS_PER_SEC;
                for(int k=0; k<testLoop; k++){
                        //呼叫 kernel (此為單一 block 的版本)
                        bubble<<<1,BLOCK>>>(gc,ga);

                        //同步化執行緒, 避免還沒做完, 量到不正確的時間
                        cudaThreadSynchronize();
                }
        t0=((double)clock()/CLOCKS_PER_SEC-t0)/testLoop;

        //測試 host 效能
        double t1=(double)clock()/CLOCKS_PER_SEC;
                for(int k=0; k<testLoop; k++){
                        bubble_host(b,a);
                }
        t1=((double)clock()/CLOCKS_PER_SEC-t1)/testLoop;

        //顯示計算時間, 並比較
        printf("time[gpu]: %g ms\n",t0*1000);
        printf("time[host]: %g ms\n",t1*1000);
        printf("ratio: %g x\n",t1/t0);

        //讀出 device 資料
        cudaMemcpy(c, gc, SIZE, cudaMemcpyDeviceToHost);


        //測試 device 結果的正確性
        printf("------------------------\n");
        bool flag=true;
        for(int k=0; k<N; k++){
                if(c[k]!=k){
                        flag=false;
                        break;
                }
        }
        printf("test[gpu]: %s\n",flag?"pass":"fail");

        //測試 host 結果的正確性
        flag=true;
        for(int k=0; k<N; k++){
                if(b[k]!=k){
                        flag=false;
                        break;
                }
        }
        printf("test[host]: %s\n",flag?"pass":"fail");

        //釋放記憶體
        cudaFree(ga);
        cudaFree(gc);
        free(a);
        free(b);
        free(c);


        return 0;
}
========================================================================

--
※ 發信站: 批踢踢實業坊(ptt.cc) 
◆ From: 114.45.213.247
推 lavatar:頭推!!
※ a5000ml:轉錄至看板 C_and_CPP                                    10/29 21:05
推 Dissipate:推@@
推 VictorTom:推
推 Luciferspear:四推
推 wch6858:我推
※ uf2000uf:轉錄至某隱形看板                                       10/30 09:58
推 moonshaped:推推推
推 vip82:感謝分享!
推 tngduh:  推