作者a5000ml (咖啡里的海洋蓝)
看板VideoCard
标题[分享] CUDA 程式设计(5) -- 第一支程式 (向量加法)
时间Wed Oct 15 20:47:02 2008
※ 第五章 第一支程式 (向量加法)
这单元以向量加法为例, 示范 CUDA 的写作, 总共有四个版本, 使用随机向量,
测试准确度和效能, 其中准确度我们拿另一个 host 的版本计算结果做为比较,
效能我们测量 kernel 和 host 函数执行数次的平均时间来比较, 分成以下四个范例
(ps. 副档名记得要 .cu 而非 .cpp)
(1) 单一区块, 单一执行绪
可以用来测试单一串流处理器 (SP, stream processor) 的效能, 并让我们容易
和传统 C 接轨, 因为只使用单一处理器, 程式码完全相同.
(2) 单一区块, 多执行绪
可以用来测试单一多处理器 (MP, multiprocessor) 的效能, 一般而言单一区块
不论包含多少个执行绪, 都使用一个 MP 来执行, 也就是 MP 里的 8 个 SP
轮流在这些执行绪间切换, 硬体上以 warp (32 个执行绪)为单位进行群组切换,
当有多个区块存在且资源足够时(暂存器&共享记忆体), 可能会有多个区块合并
到同一个 MP 执行.
(3) 多区块, 多执行绪 (不使用回圈, 用网格与区块设定代替)
当回圈数目小於 max(gridDim)*max(blockDim) 而且前後无相依性时, 可以
将回圈打平, 使用网格与区块设定代替, 每一个执行绪只计算一次, 整个核心
就是一个平行化的大回圈, 从 compute 1.0 到 1.3, 这些数目没什麽更动
每个网格最大的区块数量 max(gridDim) 是 65,535
每个区块最大的执行绪数量 max(blockDim)是 512
这种做法是为了避免回圈拖慢效能, 因为 GPGPU 的分支指令(branch, 条件判断)
是比较弱的, 理由是其执行上是以 warp 为单位, 当 32 个执行绪条件不同时,
就会导致部份的执行绪停摆, 等到下一执行周期再执行, 这种变成循序执行的
动作称为发散 (divergence), 会造成这段指令执行需要两倍时间而拖慢效能.
这种去除条件判断和回圈的程式设计法我称为「乾式(dry)」的设计法, 因为
缺乏流程控制 (flow control), 其缺点为在比较复杂的程式中容易失去弹性,
而且必需付出计算资料位址的额外成本 (每个执行绪都必需计算一次).
(4) 多区块, 多执行绪 (使用回圈)
一般而言, 回圈会造成条件发散的只有最後一个, 也就是大部份的 warp 执行
仍是良好的, 所以我们也不用那麽紧张的将回圈弄成乾式, 即使存在时常发生
发散的条件判断也不用担心, 只要把条件所附带的程式码尽量简短就行了 (让
发散变成两倍的程式码执行时间缩短, 额外的时间便能缩短), 融合条件判断
与平行化可使程式兼具弹性与效能。
===========================================================
test1.cu (单一区块, 单一执行绪)
===========================================================
以下是范例(1)的程式, 只使用到单一处理器, 核心 gpu_add() 和对照组的 host_add()
完全相同, 将它存成 test1.cu 後, 在 linux 下使用 nvcc test1.cu -O3 即可编译,
Windows 上可以修改 sdk 的 template project 来编译, 其行为如下图
执行绪 (一个一个执行)
(0) --> 记忆体
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
在 GTX 280 上执行的结果和 Intel Q9400 的比较
vector add N(1048576) elements, diff = 0
host time: 3.3 ms
gpu time: 505.8 ms
也就是单一核心 GPU 慢上 150 倍, 原因除了单核 GPU 本身计算力较弱外, 主要还是
单核心在记忆体 I/O 上的问题 (将在後续章节再讨论).
===========================================================
test1.cu 程式码
===========================================================
#include <cuda.h>
#include <time.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
//----------------------------------------------
//向量加法的运算核心 (GPU) **函式前加 __global__ 即为核心, 核心只传回 void**
__global__ void gpu_add(float* c, float* a, float* b, int n){
for(int k=0; k<n; k++){
c[k]=a[k]+b[k];
}
}
//----------------------------------------------
//向量加法的一般函式 (Host)
void host_add(float* c, float* a, float* b, int n){
for(int k=0; k<n; k++){
c[k]=a[k]+b[k];
}
}
//----------------------------------------------
//计算误差用的函式
double diff(float* a, float* b, int n){
double s=0, r=0;
for(int k=0; k<n; k++){
double w=a[k]-b[k];
s+=w*w;
r+=a[k]*a[k];
}
return sqrt(s/r); //相对误差
}
//----------------------------------------------
//时间函数 (传回单位:千分之一秒)
double ms_time(){
return (double)clock()/CLOCKS_PER_SEC*1000.0;
}
//----------------------------------------------
//主程式
int main(){
//设定向量大小
int n=1024*1024;
int size=n*sizeof(float);
//网格与区块设定
int grid=1; //gridDim (每个网格具有的区块数)
int block=1; //blockDim (每个区块具有的执行绪数)
//设定呼叫次数 (测量平均效能)
int loop=100;
//配置主机记忆体
float *a,*b,*c,*d;
a=(float*)malloc(size);
b=(float*)malloc(size);
c=(float*)malloc(size);
d=(float*)malloc(size);
//设定乱数的输入向量
srand(time(0));
for(int k=0; k<n; k++){
a[k]=(float)rand()/RAND_MAX*2-1;
b[k]=(float)rand()/RAND_MAX*2-1;
}
//配置显示卡记忆体
float *ga,*gb,*gc;
cudaMalloc((void**)&ga, size);
cudaMalloc((void**)&gb, size);
cudaMalloc((void**)&gc, size);
//载入向量 a,b 到显示卡记忆体中
cudaMemcpy(ga, a, size, cudaMemcpyHostToDevice);
cudaMemcpy(gb, b, size, cudaMemcpyHostToDevice);
//---- part 1 : 测量精确度 --------
//呼叫核心来运算 (GPU)
gpu_add<<<grid, block>>>(gc, ga, gb, n);
//呼叫一般函数来运算 (Host)
host_add(d, a, b, n);
//把计算结果存回主机
cudaMemcpy(c, gc, size, cudaMemcpyDeviceToHost);
//比较两者差异
printf("vector add N(%d) elements, diff = %g\n", n, diff(c,d,n));
//---- part 2 : 测量效能 --------
//测量 GPU 核心效能
double gpu_dt = ms_time();
for(int w=0; w<loop; w++){
gpu_add<<<grid, block>>>(gc, ga, gb, n);
cudaThreadSynchronize(); //避免核心执行不完全
}
gpu_dt = (ms_time()-gpu_dt)/loop; //平均时间
//测量 Host 函数效能
double host_dt = ms_time();
for(int w=0; w<loop; w++){
host_add(d, a, b, n);
}
host_dt = (ms_time()-host_dt)/loop; //平均时间
//输出平均执行时间
printf("host time: %g ms\n", host_dt);
printf("gpu time: %g ms\n", gpu_dt);
//释放主机记忆体
free(a);
free(b);
free(c);
free(d);
//释放显示卡记忆体
cudaFree(ga);
cudaFree(gb);
cudaFree(gc);
return 0;
}
===========================================================
test2.cu (单一区块, 多执行绪)
===========================================================
以下是范例(2)的需要修改的部份, 其中我们用到两个内建唯读变数
threadIdx.x 区块中的执行绪索引
blockDim.x 区块中的执行绪数量
假设我们有 8 个执行绪, 其 id 为 0~7, 这些平行的执行绪可以想向成一排印章,
沿着记忆体一块一块的印过去, 每次步进 8 个单位, 如下图所示
执行绪区块
(01234567) --> 记忆体
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
(01234567) 步进 8 个单位的距离 (blockDim.x)
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
(01234567) 再步进 8 个单位
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
...
(01234567)
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
执行的结果为
vector add N(1048576) elements, diff = 0
host time: 3.7 ms
gpu time: 12.2 ms
还是比 CPU 慢上一些, 因为还有很多 MP 没有发挥.
===========================================================
test2.cu 程式码在 test1.cu 要修改的部份
===========================================================
//运算核心
__global__ void gpu_add(float* c, float* a, float* b, int n){
for(int k=threadIdx.x; k<n; k+=blockDim.x){
c[k]=a[k]+b[k];
}
}
int main(){
...
//网格与区块设定
int grid=1; //只使用一个区块
int block=512;
...
}
===========================================================
test3.cu (多区块, 多执行绪) 不使用回圈,用网格与区块设定代替
===========================================================
再来是使用乾式设计法的范例(3), 先定出 global idx, 然後一次印完整个向量,
使用到三个内建唯读变数
blockIdx.x 网格中的区块 id
blockDim.x 区块中的执行绪数量
threadIdx.x 区块中的执行绪 id
其中全域索引可用下式定址
global idx = blockIdx.x*blockDim.x + threadIdx.x
行为如下
global idx
(012345678...............................................N)
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
记忆体
执行结果为
vector add N(1048576) elements, diff = 0
host time: 3.56 ms
gpu time: 0.12 ms
所以 GPU 的能力被释放出来了, 比 CPU 快上 30 倍.
===========================================================
test3.cu 程式码在 test1.cu 要修改的部份
===========================================================
__global__ void gpu_add(float* c, float* a, float* b, int n){
int j=blockIdx.x*blockDim.x+threadIdx.x;
c[j]=a[j]+b[j];
}
int main(){
...
//网格与区块设定
int block=512;
int grid=n/block; //必需能被整除
...
}
===========================================================
test4.cu (多区块, 多执行绪) 使用回圈,较为一般性
===========================================================
再来是使用回圈的范例(4), 先定出 global id, 但是它并不像范例(3)一样包含整个
向量大小, 然後像范例(2) 一块一块印过去
行为如下
global id (M = gridDim*blockDim)
(01234567...M-1) --> 记忆体
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
(01234567...M-1) 步进 M 个单位
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
...
(01234567...M-1)
ooooooooooooooooooooooooooooooooooooooooooooooooooooooooo
执行结果为
vector add N(1048576) elements, diff = 0
host time: 3.64 ms
gpu time: 0.12 ms
和乾式的效能差不多, 但向量长度变得更有弹性, 可以和 blockDim & gridDim 无关.
===========================================================
test4.cu 程式码在 test1.cu 要修改的部份
===========================================================
__global__ void gpu_add(float* c, float* a, float* b, int n){
int j=blockIdx.x*blockDim.x+threadIdx.x;
int m=gridDim.x*blockDim.x;
for(int k=j; k<n; k+=m){
c[k]=a[k]+b[k];
}
}
int main(){
...
//网格与区块设定
int block=256;
int grid=30;
...
}
===========================================================
※ 使用到的 CUDA API ※
---------------------------------------------------------------
(1)显示卡记忆体配置 cudaMalloc()
cudaError_t cudaMalloc(void** ptr, size_t count);
ptr 指向目的指位器之位址
count 欲配置的大小(单位 bytes)
传回值 cudaError_t 是个 enum, 执行成功时传回 0, 其它的错误代号可用
cudaGetErrorString() 来解译.
---------------------------------------------------------------
(2)显示卡记忆体释放 cudaFree()
cudaError_t cudaFree(void* ptr);
ptr 指向欲释放的位址 (device memory)
---------------------------------------------------------------
(3)记忆体拷贝 cudaMemcpy()
cudaError_t cudaMemcpy(void* dst, const void* src, size_t count,
enum cudaMemcpyKind kind);
dst 指向目的位址
src 指向来源位址
count 拷贝区块大小 (单位 bytes)
kind 有四种拷贝流向
cudaMemcpyHostToHost 主机 -> 主机
cudaMemcpyHostToDevice 主机 -> 装置
cudaMemcpyDeviceToHost 装置 -> 主机
cudaMemcpyDeviceToDevice 装置 -> 装置
---------------------------------------------------------------
(4) 错误字串解译 cudaGetErrorString()
const char* cudaGetErrorString(cudaError_t error);
传回错误代号(error)所代表的字串
---------------------------------------------------------------
※名词解释※
(1) 执行缆 (warp): 包含 32 个执行绪的单元, 多处理器(MP) 的执行绪切换以此为单位.
(2) 乾式 (dry): 相对於 flow 而言, 缺乏流程控制 (flow control),使用上缺乏弹性.
(3) 发散 (divergence): 同一 warp 的执行绪条件不同时, 就会导致部份的执行绪停摆,
等到下一执行周期再循序执行.
(4) 全域索引 (global idx): 由 (blockIdx,threadIdx) 合成的整体索引.
--
※ 发信站: 批踢踢实业坊(ptt.cc)
◆ From: 114.45.209.219
※ a5000ml:转录至看板 C_and_CPP 10/15 21:00
※ 编辑: a5000ml 来自: 114.45.209.219 (10/15 21:07)
1F:推 dkfum:快M 10/15 21:28
2F:推 dsin:请问一下 这有要什麽编译器(VC6 VC2008) 和什麽lib才能用吗? 10/15 22:00
3F:推 wowtiger:.NET 2002 .NET 2003 限定 10/15 22:00
4F:→ dsin:唷唷 翻旧文翻到了 抱歉 10/15 22:01
5F:推 ycjcsie:我用2005也OK耶 10/15 22:08
6F:推 kuninaka:有LINUX版本吗@@ 10/15 22:09
7F:推 sdk:回楼上 有 10/15 22:30
8F:推 Smener:好文推~~ 10/15 23:30
9F:推 Dissipate:专业 10/16 00:28
10F:推 Luciferspear:9 10/16 00:45
11F:推 mnmnqq:这系列的文都有收喔~ 10/16 00:53
12F:→ GelionLin:明灯啊!!! 10/16 01:32
13F:推 f7258:快推!! 以免大家说我看不懂XD 10/16 09:07
14F:推 VictorTom:推.... 10/16 10:31
15F:推 lookers:同楼楼上XD 10/16 10:34
16F:推 moonshaped:好东西要推 要继续下去~~~~ 10/16 14:28
17F:推 vixen:太专业了 10/16 14:37
18F:推 lavatar:好文当推!! 10/16 15:01
19F:推 sardine:看两页就打哈欠了 10/16 19:33
20F:推 xxxEVA:终於开始看不懂了XD 10/16 19:47
※ uf2000uf:转录至某隐形看板 10/16 21:49
21F:推 finkel:没记错的话,warp内遇到divergence branch时,其实if、else 10/18 02:21
22F:→ finkel:都会去执行的,只是最後会只留正确的branch 10/18 02:22
23F:→ a5000ml:是啊, 只是执行时不合条件的 thread 会被 disable 10/23 06:12
24F:推 finkel:我会降说是有想过,若都会执行,那在某些时候写else if不就 10/25 01:18
25F:→ finkel:没有实质意义 10/25 01:18
26F:→ a5000ml:它也不是都会执行, 没 diverge 就只执行一个 branch, 10/25 11:09
27F:→ a5000ml:另一个是空的, 若 diverge 时, 两种都会执行, 不合条件的 10/25 11:10
28F:→ a5000ml:threads 会被 disable, SP 会略过它, 其实它的 branch 和 10/25 11:12
29F:→ a5000ml:单核心的一样, 只是变批次处理而己, 要执行的 threads数量 10/25 11:15
30F:→ a5000ml:比较少的 branch 还是会快一些 10/25 11:16