CUDA: GPU高性能运算

时间：2021-08-04 15:49:52

CUDA: GPU高性能运算

 2013-10-11 22:23 5650人阅读 评论(0)收藏举报 CUDA: GPU高性能运算分类：

目录(?)[+]

0 序言

CUDA是异构编程的一个大头，洋洋洒洒的看了些资料，但是，感觉这个技术没有像C++或者Java那样有自己的权威的《编程思想》来指导系统学习，总是感觉心里不踏实，是不是自己还没掌握深入、或者说心里没底气说自己已经入门了、已经熟悉了、已经精通了。站在一个初学者的角度，作为一个笔记式的记录，讲解自己学习和理解CUDA过程中的一些列想到的、碰到的问题。享受一个东西不一定是结果，可以是从无知到了解到精通的这个整个过程。

1 给自己提几个问题

对的，我想要做什么事情的时候，习惯性的给自己提如下问题：

问题1：GPU高性能运算之CUDA（下午简称CUDA技术）是干嘛用的，我为什么要学它？

这个问题如果我回答是，纯粹为了掌握一门新的技术，如果你还是学生则可以，如果你是一个手里有项目的工程师，那么，我觉得没什么必要去学这个东西。我个人理解CUDA是计算机里面的边缘技术，是对程序执行性能的提高的一种方式，可以理解成一个工具。工具这种东西，你需要用的时候再去学，你不需要用的时候，你知道有这回事情就可以了。如果你用不到的东西，你也很贪心什么都去学，第一学不精，第二计算机的技术太多太广，对应我这样的智商一般的人来说是不靠谱的。

我要学习CUDA，因为项目的任务可分解性较大，粒度相关性小，某项目大概一个计算任务可以分解成3000*3000*300规模的计算，而且这些元素之间完全独立。又加入，你需要对全市100万人，每个人计算下该个体的每年的收入与支出的净值（从1900年计算到2000年），淡然有人说放excel中计算不就得了，好吧……我只是讲个通俗点的例子。那么这个计算规模是1000000*10，倘若要对每个人每年在做点其他什么高级的算法得出一个什么指数，恐怕excel还是实现起来不太容易。所以，这种并行度很大的问题，我们可以考虑用CUDA来解决。也许你也知道，搞图像的，CUDA就显得很重要啦。


问题2：我的基础是什么？

学习新的技术，我喜欢和之前学过的某个还熟悉的东西做为比对，这样理解起来可能会快一些。比如学Java的时候，我想着C++，学UML的时候，我想着面向对象编程，学CUDA呢？我有什么么？可能很多人会没有异构编程的历史。我很幸运，之前弄过半年的OpenMP，对的——多核编程技术。后来我了解到，如果是多个GPU，是需要用到OpenMP来做的。OpenMP本身和CUDA好像并没什么关系，但是，片上多核的并行算法是很想通的。很好的一个例子是：奇偶排序。哈哈，给自己提了点学好CUDA的信心。

好了，明白自己的需求和基础，给自己一个学习的定位，什么地方该花时间去琢磨，相比应该很清楚。

2 你好，GPU！

GPU有两大开发商——英伟达和AMD，支持CUDA的是英伟达，很好啊，之前买电脑是有先见之明的，买了英伟达的显卡——GT520，不是特别高端，和单位的GTX650ti2G比起来，有那么一点逊色，但是，好歹可以跑CUDA！

英伟达支持CUDA编程的显卡型号从G8800开始，都是可以的。一开始作为图像处理用，而今，天文地理、数学金融、医疗军事等等，都开始尝试发挥GPU的优势。

GPU的计算核心也是隔年换代，现在已经倒GK10X了，计算能力也是逐渐提升，目前已经最好的有3.0。GPU的架构从原始的到费米的，再到开普勒的，我们没必要去一个个了解，我们先了解GPU的这个大概的历史，免得和人交谈说不出一和二，脱离菜鸟嘛！

GPU和CPU的区别可以参考下，讲的还算详细和明了。

http://www.cnblogs.com/viviman/archive/2012/11/26/2789113.html

可以这样的去理解，单核CPU多线程并行，是感官上的并行，世界上在CPU上还是串行指令在跑；而在GPU上，才叫真正的并行！需要介绍下CUDA1.0开发包是支持在CPU上模拟GPU开发的，其原理就是用多线程来模拟；而现在的版本就不支持了。最新的是5.0的，官网上是下不到之前的了，反正我是找到腿软了还没找到。1.0是古董了！如果你有，一定要给我开开眼见哦。

我们的CUDA编程，很明显是GPU和CPU一起来处理的嘛——异构编程！对的，我自己的理解是就一个工程而言：CPU处理串行计算业务，GPU处理并行计算业务，这里将的并行都是并行度可观的哦，不是说两个元素你也来GPU上计算，那样是不环保的——会浪费很多GPU的资源！

说道环保，我想多说一句，在GPU编程，很体现“绿色”理念。48个核，你如果写的好，48个核全在干活，而且干的是有意义的活，那么你是合理的利用资源，如果你只让一个核在干有意义的活，其他的都在空转，那么你很浪费电哦。

对于CPU和GPU分配自己的业务，稍微画个图失意一下，如图1：

CUDA: GPU高性能运算

图1 工程中GPU和CPU的分工

总的来说，GPU只是干计算并行度高的功能模块的活，一定不可以越权啦！

3 你好，CUDA！

3.1 开发环境配置

第一次和同事交谈，我说你和我说说C-U-D-A。他说：哭打……苦打…… ……。半天后，我说，你和我说说C-U-D-A，你说哭打是什么东西……它说哭打就是C-U-D-A。我操，顿时傻眼了，原来这东西行业里年哭打，对的，Cu - Da，连起来就是这么发音的，我的无知啊。我们还是用中文解释吧：CUDA的意思是统一计算设备架构。

CUDA的集成开发环境可以参考下：

http://www.cnblogs.com/viviman/archive/2012/11/05/2775100.html

在win7+vs2008+CUDA5.0的环境下体验，是一种新的尝试，我自己配的时候，很少或者就没有5.0的配置博客文档等等。因为5.0的那一场雪比2010年来的更晚一些……

5.0是和2.0、3.0、4.0都会有那么一点不同的，是集成了SDK和TOOL两个东西，之前是分开的，现在是合在一起的。其实是差不多的，但是，这一分一合，就会给人很不习惯。不过，我虽然愚钝，但是，试了几次之后还是摸索成功。当第一个helloworld输出后，心里是有那么一点小激动的，立马跑到小区门口，买了半斤羊肉吃了，因为为了配这一套环境，我中午饭都没吃！我这只哭打小菜鸟就是可怜啊！

3.2 特殊的"hello world"

搭建好环境，你肯定想看看我的helloworld程序，对的，但是我不想输出helloworld，我想干一件事情是：我在CPU上创建个变量，传到GPU中，然后，在GPU中赋值，然后传出来。这件事情如果成功了，是不是可以说明，通了+GPU工作了！网上一搜的CUDA的helloworld程序，都是在CPU上输出helloworld，那多不过瘾。

[cpp] view plain copy
  1. __global__ void hello(char *ch)  
  2. {  
  3. ch = {'h', 'e', 'l', 'l', 'o'};  
  4. }  
  5. int main()  
  6. {  
  7. ……  
  8. hello<<<1,1>>>(dev_ch);  
  9. ……  
  10. return 0;  
  11. }  

I think you know my idea.

在你网上搜索到的程序的基础上，作这样的一个改变，相信自己动手的才是快乐的。

4 敲开编程的门

我习惯性的喜欢先看一门语言的关键字，CUDA的关键字很简单很少：

函数类型

__global__: 用来修饰内核函数的，内核函数是什么呢，内核函数是跑在GPU上的函数；与之对应的是主机函数，用__host__修饰，也可以缺省，跑在CPU上。因此，CPU也叫主机，GPU也叫设备。通常定义这个内核函数，我喜欢在函数名前加个kernel作为修饰，让自己清楚点。

比如__global__ void kerneladd(float *a){}

__device__: 也是用来修饰内核函数的，那和__global__有什么区别吗？对的。__global__修饰的内核函数只能被主机函数调用；__device__修饰的内核函数只能被内核函数调用，应该很好理解。

__host__: 主机函数，供主机函数调用，可缺省哦，一般情况下，都是缺省的，知道这个东西就行。

存储类型

寄存器：在核函数内 int i即表示寄存器变量。

__global__：     全局内存。在主机函数中开辟和释放。

__shared__：   共享存储，每个block内的线程共享这个存储。

__constant__：常量存储，只读。定义在所有函数之外，作用范围整个文件。

__texture__：   纹理存储，只读。内存不连续。

CUDA: GPU高性能运算

CUDA: GPU高性能运算

内建变量

dim3, threadId, blockId,  gridId

5 GPU也不允许偏心

并行的事情多了，我们作为GPU的指令分配者，不能偏心了——给甲做的事情多，而乙没事做，个么甲肯定不爽的来。所以，在GPU中，叫做线程网络的分配。首先还是来看下GPU的线程网络吧，图2：

CUDA: GPU高性能运算

图2 线程网络

我们将具体点的，在主机函数中如果我们分配的是这样的一个东西：

dim3 blocks(32, 32);

dim3 threads(16, 16);

dim3是神马？dim3是一个内置的结构体，和linux下定义的线程结构体是个类似的意义的东西，dim3结构变量有x，y，z，表示3维的维度。不理解没关系，慢慢看。

kernelfun<<<blocks, threads>>>();

我们调用kernelfun这个内核函数，将blocks和threads传到<<<,>>>里去，这句话可牛逼大了——相当于发号施令，命令那些线程去干活。这里使用了32*32 * 16*16个线程来干活。你看明白了吗？blocks表示用了二维的32*32个block组，而每个block中又用了16*16的二维的thread组。好吧，我们这个施令动用了262144个线程！我们先不管GPU内部是如何调度这些线程的，反正我们这一句话就是用了这么多线程。

那我们的内核函数kernelfun()如何知道自己执行的是哪个线程？这就是线程网络的特点啦，为什么叫网络，是有讲究的，网络就可以定格到网点：

比如int tid = threadId.x + blockId.x * 16

这里有一个讲究，block是有维度的，一维、二维、三维。

对于一维的block，tid = threadId.x

对于（Dx，Dy）二维的block，tid = threadId.x + Dx*threadId.y

对于（Dx，Dy，Dz）三维的block，tid = threadId.x + Dx*threadId.y + Dz*Dy*threadId.z

我习惯的用这样的模式去分配，比较通用：

dim3 dimGrid();

dim3 dimBlock();

kerneladd<<<dimGrid, dimBlock>>>();

这可是万金油啊，你需要做的事情是填充dimGrid和dimBlock的结构体构造函数变量，比如，dimGrid(16, 16)表示用了16*16的二维的block线程块。

（0,0）（0,1）（0,2）……（0,15）

（1,0）（1,1）（1,2）……（1,15）

（2,0）（2,1）（2,2）……（2,15）

……

（15,0）（15,1）（15,2）……（15,15）

（，）是（dimGrid.x， dimGrid.y）的网格编号。

我们这么理解吧，现在又一群人，我们分成16*16个小组（block），排列好，比如第3行第4列就指的是（2,3）这个小组。

而dimBlock（16,16）表示每个小组有16*16个成员，如果你想点名第3行第4列这个小组的里面的第3行第4列那个同学，那么，你就是在（2,3）这个block中选择了（2,3）这个线程。这样应该有那么一点可以理解进去的意思了吧？不理解透彻么什么关系，这个东西本来就是cuda中最让我纠结的事情。我们且不管如何分配线程，能达到最优化，我们的目标是先让GPU正确地跑起来，计算出结果即可，管他高效不高效，管他环保不环保。

唠叨了这么多，下面我们用一个最能说明问题的例子来进一步理解线程网络分配机制来了解线程网络的使用。

一维网络线程

eg：int arr[1000]，对每个数组元素进行加1操作。

idea：我们最直接的想法，是调度1000个线程去干这件事情。

first pro：我想用一个小组的1000个人员去干活。这里会存在这样一个问题——一个小组是不是有这么多人员呢？是的，这个事情你必须了解，连自己组内多少人都不知道，你也不配作指挥官呀。对的，这个参数叫做maxThreadsPerBlock，如何取得呢？

好吧，cuda定义了一个结构体cudaDeviceProp，里面存入了一系列的结构体变量作为GPU的参数，出了maxThreadsPerBlock，还有很多信息哦，我们用到了再说。

maxThreadsPerBlock这个参数值是随着GPU级别有递增的，早起的显卡可能512个线程，我的GT520可以跑1024个线程，办公室的GTX650ti2G可以跑1536个，无可非议，当然多多益善。一开始，我在想，是不是程序将每个block开的线程开满是最好的呢？这个问题留在以后在说，一口吃不成胖子啦。

好吧，我们的数组元素1000个，是可以在一个block中干完的。

内核函数：

[cpp] view plain copy
  1. #define N 1000  
  2. __gloabl__ void kerneladd(int *dev_arr)  
  3. {  
  4.     int tid = threadId.x;  
  5.     if (tid < 1000)  
  6.     dev_arr[tid] ++;  
  7. }  
  8. int main()  
  9. {  
  10.     int *arr, *dev_arr;// 习惯的我喜欢在内核函数参数变量前加个dev_作为标示  
  11.     // 开辟主机内存，arr = (int*)malloc(N*sizeof(int));  
  12.     // 开辟设备内存  
  13.     // 主机拷贝到设备  
  14.     kerneladd<<<1, N>>>(dev_arr);  
  15.     // 设备拷贝到主机  
  16.     // 打印  
  17.     // 释放设备内存  
  18.     // 释放主机内存  
  19.     return 0;  
  20. }  

呀，原来这么简单，个么CUDA也忒简单了哇！这中想法是好的，给自己提高信心，但是这种想法多了是不好的，因为后面的问题多了去了。

盆友说，1000个元素，还不如CPU来的快，对的，很多情况下，数据量并行度不是特别大的情况下，可能CPU来的更快一些，比较设备与主机之间互相调度操作，是会有额外开销的。有人就问了，一个10000个元素的数组是不是上面提供的idea就解决不了啦？对，一个block人都没怎么多，如何完成！这个情况下有两条路可以选择——

第一，我就用一个组的1000人来干活话，每个人让他干10个元素好了。

这个解决方案，我们需要修改的是内核函数：

[cpp] view plain copy
  1. __global__ void kernelarr(int *dev_arr)  
  2. {  
  3.     int tid = threadId.x;  
  4.     if(tid < 1000) // 只用0~999号线程  
  5.     {  
  6.         //每个线程处理10个元素，比如0号线程处理0、1001、2001、……9001  
  7.         for(int i = tid; i<N; i=i+1000)  
  8.         {  
  9.         dev_arr[tid] ++;  
  10.         }  
  11.     }  
  12. }  

第二，我多用几个组来干这件事情，比如我用10个组，每个组用1000人。

这个解决方案就稍微复杂了一点，注意只是一点点哦~因为，组内部怎么干活和最原始的做法是一样的，不同之处是，我们调遣了10个组去干这件事情。

首先我们来修改我们的主机函数：

[cpp] view plain copy
  1. int main()  
  2. {  
  3. ……  
  4. kerneladd<<<10, 1000>>>(dev_arr);//我们调遣了10个组，每个组用了1000人  
  5. ……  
  6. }  

盆友要问了，10个组每个组1000人，你怎么点兵呢？很简单啊，第1组第3个线程出列，第9组第9个线程出列。每个人用组号和组内的编号定了位置。在线程网络中，blockId.x和threadId.x就是对应的组号和组内编号啦，我必须要这里开始形象点表示这个对应关系，如果这个对应关系是这样子的[blockId.x，threadId.x]，那么我们的数组arr[10000]可以这样分配给这10个组去干活：

(0,0)——arr[0]，(0,1)——arr[1]，……(0,999)——arr[999]

(1,0)——arr[0+1*1000]，(1,1)——arr[1+1*1000]，……(1,999)——arr[999+1*1000]

……

(9,0)——arr[0+9*1000]，(9,1)——arr[1+9*1000]，……(9,999)——arr[999+9*1000]

是不是很有规律呢？对的，用blockId.x和threadId.x可以很好的知道哪个线程干哪个元素，这个元素的下表就是threadId.x + 1000*blockId.x。

这里我想说的是，如果我们哪天糊涂了，画一画这个对应关系的表，也许，就更加清楚的知道我们分配的线程对应的处理那些东西啦。 

一维线程网络，就先学这么多了。

二维网络线程

eg2：int arr[32][16]二维的数组自增1。

第一个念头，开个32*16个线程好了哇，万事大吉！好吧。但是，朕现在想用二维线程网络来解决，因为朕觉得一个二维的网络去映射一个二维的数组，朕看的更加明了，看不清楚自己的士兵，如何带兵打仗！

我还是画个映射关系：

一个block中，现在是一个二维的thread网络，如果我用了16*16个线程。

(0,0)，(0,1)，……(0,15)

(1,0)，(1,1)，……(1,15)

……

(15,0)，(15,1)，……(15,15)

呀，现在一个组内的人称呼变了嘛，一维网络中，你走到一个小组里，叫3号出列，就出来一个，你现在只是叫3号，没人会出来！这个场景是这样的，现在你班上有两个人同名的人，你只叫名，他们不知道叫谁，你必须叫完整点，把他们的姓也叫出来。所以，二维网络中的(0,3)就是原来一维网络中的3，二维中的(i,j)就是一维中的(j+i*16)。不管怎么样，一个block里面能处理的线程数量总和还是不变的。

一个grid中，block也可以是二维的，一个block中已经用了16*16的thread了，那我们一共就32*16个元素，我们用2个block就行了。

先给出一个代码清单吧，程序员都喜欢看代码，这段代码是我抄袭的。第一次这么完整的放上代码，因为我觉得这个代码可以让我说明我想说的几个问题：

第一，二维数组和二维指针的联系。

第二，二维线程网络。

第三，cuda的一些内存操作，和返回值的判断。

[cpp] view plain copy
  1. #include <stdio.h>   
  2. #include <stdlib.h>   
  3. #include <cuda_runtime.h>   
  4.    
  5. #define ROWS 32   
  6. #define COLS 16   
  7. #define CHECK(res) if(res!=cudaSuccess){exit(-1);}   
  8. __global__ void Kerneltest(int **da, unsigned int rows, unsigned int cols)   
  9. {   
  10.  unsigned int row = blockDim.y * blockIdx.y + threadIdx.y;   
  11.  unsigned int col = blockDim.x * blockIdx.x + threadIdx.x;   
  12.  if (row < rows && col < cols)   
  13.  {   
  14.   da[row][col] = row * cols + col;   
  15.  }   
  16.  }   
  17.    
  18. int main(int argc, char **argv)   
  19. {   
  20.  int **da = NULL;   
  21.  int **ha = NULL;   
  22.  int *dc = NULL;   
  23.  int *hc = NULL;   
  24.  cudaError_t res;   
  25.  int r, c;   
  26.  bool is_right = true;   
  27.    
  28.  res = cudaMalloc((void**)(&da), ROWS * sizeof(int*)); CHECK(res)   
  29.  res = cudaMalloc((void**)(&dc), ROWS * COLS * sizeof(int)); CHECK(res)   
  30.  ha = (int**)malloc(ROWS * sizeof(int*));   
  31.  hc = (int*)malloc(ROWS * COLS * sizeof(int));   
  32.    
  33.  for (r = 0; r < ROWS; r++)   
  34.  {   
  35.   ha[r] = dc + r * COLS;   
  36.  }   
  37.  res = cudaMemcpy((void*)(da), (void*)(ha), ROWS * sizeof(int*), cudaMemcpyHostToDevice); CHECK(res)   
  38.  dim3 dimBlock(16, 16);   
  39.  dim3 dimGrid((COLS + dimBlock.x - 1) / (dimBlock.x), (ROWS + dimBlock.y - 1) / (dimBlock.y));   
  40.  Kerneltest<<<dimGrid, dimBlock>>>(da, ROWS, COLS);   
  41.  res = cudaMemcpy((void*)(hc), (void*)(dc), ROWS * COLS * sizeof(int), cudaMemcpyDeviceToHost); CHECK(res)   
  42.    
  43.  for (r = 0; r < ROWS; r++)   
  44.  {   
  45.   for (c = 0; c < COLS; c++)   
  46.   {   
  47.    printf("%4d ", hc[r * COLS + c]);   
  48.    if (hc[r * COLS + c] != (r * COLS + c))   
  49.    {   
  50.     is_right = false;   
  51.    }   
  52.   }   
  53.   printf("\n");   
  54.  }   
  55.  printf("the result is %s!\n", is_right ? "right" : "false");   
  56.  cudaFree((void*)da);   
  57.  cudaFree((void*)dc);   
  58.  free(ha);   
  59.  free(hc);   
  60.  getchar();   
  61.  return 0;   
  62.  }   

简要的来学习一下二维网络这个知识点， 

dim3 dimBlock(16, 16); 

//定义block内的thread二维网络为16*16

dim3 dimGrid((COLS + dimBlock.x - 1) / (dimBlock.x),  (ROWS + dimBlock.y - 1) / (dimBlock.y)); 

//定义grid内的block二维网络为1*2

unsigned int row = blockDim.y * blockIdx.y  +  threadIdx.y; 

//二维数组中的行号

unsigned int col = blockDim.x * blockIdx.x  +  threadIdx.x; 

//二维线程中的列号

三维网络线程

dim3定义了三维的结构，但是，貌似二维之内就能处理很多事情啦，所以，我放弃学习三维。网上看到的不支持三维网络是什么意思呢？先放一放。

给自己充充电

同一块显卡，不管你是二维和三维或一维，其计算能力是固定的。比如一个block能处理1024个线程，那么，一维和二维线程网络是不是处理的线程数一样呢？

 

CUDA: GPU高性能运算

回答此问题，先给出网络配置的参数形式——<<<Dg,Db,Ns,S>>>，各个参数含义如下：

Dg：定义整个grid的维度，类型Dim3，但是实际上目前显卡支持两个维度，所以，dim3<<Dg.x, Dg.y, 1>>>第z维度默认只能为1，上面显示出这个最大有65536*65536*1，每行有65536个block，每列有65536个block，整个grid中一共有65536*65536*1个block。

Db：定义了每个block的维度，类型Dim3，比如512*512*64，这个可以定义3维尺寸，但是，这个地方是有讲究了，三个维度的积是有上限的，对于计算能力1.0、1.1的GPU，这个值不能大于768，对于1.2、1.3的不能大于1024，对于我们试一试的这块级别高点的，不能大于1536。这个值可以获取哦——maxThreadsPerBlock

Ns：这个是可选参数，设定最多能动态分配的共享内存大小，比如16k，单不需要是，这个值可以省略或写0。

S：也是可选参数，表示流号，默认为0。流这个概念我们这里不说。

接着，我想解决几个你肯定想问的两个问题，因为我看很多人想我这样的问这个问题：

1 block内的thread我们是都饱和使用吗？

答：不要，一般来说，我们开128或256个线程，二维的话就是16*16。

2 grid内一般用几个block呢？

答：牛人告诉我，一般来说是你的流处理器的4倍以上，这样效率最高。

回答这两个问题的解释，我想抄袭牛人的一段解释，解释的好的东西就要推广呀：

GPU的计算核心是以一定数量的Streaming Processor(SP)组成的处理器阵列，NV称之为Texture Processing Clusters(TPC)，每个TPC中又包含一定数量的Streaming Multi-Processor(SM)，每个SM包含8个SP。SP的主要结构为一个ALU（逻辑运算单元），一个FPU（浮点运算单元）以及一个Register File(寄存器堆)。SM内包含有一个Instruction Unit、一个Constant Memory、一个Texture Memory，8192个Register、一个16KB的Share Memory、8个Stream Processor(SP)和两个Special Function Units（SFU）。（GeForce9300M GS只拥有1个SM） Thread是CUDA模型中最基本的运行单元，执行最基本的程序指令。Block是一组协作Thread，Block内部允许共享存储，每个Block最多包含512个Thread。Grid是一组Block，共享全局内存。Kernel是在GPU上执行的核心程序，每一个Grid对应一个Kernel任务。 在程序运行的时候，实际上每32个Thread组成一个Warp，每个 warp 块都包含连续的线程，递增线程 ID 。Warp是MP的基本调度单位，每次运行的时候，由于MP数量不同，所以一个Block内的所有Thread不一定全部同时运行，但是每个Warp内的所有Thread一定同时运行。因此，我们在定义Block Size的时候应使其为Warp Size的整数倍，也就是Block Size应为32的整数倍。理论上Thread越多，就越能弥补单个Thread读取数据的latency ，但是当Thread越多，每个Thread可用的寄存器也就越少，严重的时候甚至能造成Kernel无法启动。因此每个Block最少应包含64个Thread，一般选择128或者256，具体视MP数目而定。一个MP最多可以同时运行768个Thread，但每个MP最多包含8个Block，因此要保持100%利用率，Block数目与其Size有如下几种设定方式： Ø 2 blocks x 384 threads Ø 3 blocks x 256 threads Ø 4 blocks x 192 threads Ø 6 blocks x 128 threads Ø 8 blocks x 96 threads 

这些电很重要啊，必须要充！不然，我就很难理解为什么网络线程如何分配的。

6 规约思想和同步概念

扩大点说，并行计算是有一种基本思想的，这个算法能解决很多很常规的问题，而且很实用，比如说累加和累积等——规约思想。对于基础的、重要的，我想有必要系统的学习。

我觉得有必要重新复制下之前写的这篇介绍：

http://www.cnblogs.com/viviman/archive/2012/11/21/2780286.html

并行程序的开发有其不同于单核程序的特殊性，算法是重中之重。根据不同业务设计出不同的并行算法，直接影响到程序的效率。因此，如何设计并行程序的算法，似乎成为并编程的最大难点。观其算法，包括cuda sdk的例子和网上的牛人，给出的一些例子，以矩阵和矢量处理为主，深入点的包括fft和julia等数学公式，再高级一点的算是图形处理方面的例子。学习这些算法的思想，免不了有自己的一点点总结。之前学习过omp编程，结合现在的cuda，我觉得要理解并行编程，首先理解划分和规约这两个概念。也许你的算法学的更加扎实。划分是《算法》里面的一个重要思想，将一个大的问题或任务，分解成小问题小任务，各个击破，最后归并结果；规约是《cuda**》书上介绍的一个入门的重要思想，规约算法(reduction)用来求连加、连乘、最值等，应用广泛。每次循环参加运算的线程减少一半。不管算法的思想如何花样，万变不离其中的一点--将一个大的任务分解成小的任务集合，分解原则是粒度合适尽量小、数据相关性尽量小。如此而已。因为，我们用GPU是为了加速，要加速必须提高执行任务的并行度！明白这个道理，那么我们将绞尽脑汁地去想方设法分析自己手上的任务，分解、分解、分解！这里拿规约来说事情，因为，规约这个东西，似乎可以拿来单做9*9乘法表来熟悉，熟悉了基础的口诀，那么99*99的难题也会迎刃而解。ex：矢量加法，要实现N=64*256长度的矢量的累加和。假设a+b计算一次耗时t。

cpu计算：显然单核的话需要64*256*t。我们容忍不了。

gpu计算：最初的设想，我们如果有个gpu能同时跑N/2个线程，我们这N/2个线程同时跑，那么不是只需要t时间就能将N个数相加编程N/2个数相加了吗？对的。这一轮我们用了t时间；接着的想法，我们不断的递归这个过程，能发现吗？第二轮，我们用N/2/2个线程同时跑，剩下N/2/2个数相加，这一轮我们同样用了t时间；一直这样想下去吧，最后一轮，我们用了1个线程跑，剩下1个数啦，这就是我们的结果！每一轮时间都为t，那么理想情况，我们用了多少轮这样的计算呢？计算次数=log(N)=6*8=48，对的，只用了48轮，也就是说，我们花了48*t的时间！

规约就是这样，很简单，很好用，我们且不管程序后期的优化，单从这个算法分析上来说，从时间复杂度N降到了logN，这在常规算法上，要提高成这样的效率，是不得了的，这是指数级别的效率提高！所以，你会发现，GPU有CPU无法取代的得天独厚的优势——处理单元真心多啊！

规约求和的核函数代码如下：

[cpp] view plain copy
  1. __global__ void RowSum(float* A, float* B)  
  2. {  
  3.  int bid = blockIdx.x; int tid = threadIdx.x;  
  4.   
  5.  __shared__ s_data[128]; //read data to shared memory   
  6.  s_data[tid] = A[bid*128 + tid];   
  7.  __synctheads(); //sync  
  8.  for(int i=64; i>0; i/=2)  
  9.  {  
  10.   if(tid<i) s_data[tid] = s_data[tid] + s_data[tid+i] ;  
  11.   __synctheads();  
  12.  }  
  13.  if(tid==0)  
  14.    B[bid] = s_data[0];  
  15. }  

这个例子还让我学到另一个东西——同步！我先不说同步是什么，你听我说个故事：我们调遣了10个小组从南京去日本打仗，我们的约定是，10个组可以自己行动，所有组在第三天在上海机场会合，然后一起去日本。这件事情肯定是需要处理的，不能第1组到了上海就先去日本了，这些先到的组，唯一可以做的事情是——等待！这个先来后到的事情，需要统一管理的时候，必须同步一下，在上海这个地方，大家统一下步调，快的组等等慢的组，然后一起干接下去的旅程。

是不是很好理解，这就是同步在生活中的例子，应该这样说，计算机的所有机制和算法很多都是源于生活！结合起来，理解起来会简单一点。

在CUDA中，我们的同步机制用处大吗？又是如何用的呢？我告诉你，一个正常规模的工程中，一般来说数据都会有先来后到的关系，这一个计算结果可能是提供给另一个线程用的，这种依赖关系存在，会造成同步的应用。

__synctheads()这句话的作用是，这个block中的所有线程都执行到此的时候，都听下来，等所有都执行到这个地方的时候，再往下执行。

7 撬开编程的锁

锁是数据相关性里面肯定要用到的东西，很好，生活中也一样，没锁，家里不安全；GPU中没锁，数据会被“盗”。

对于存在竞争的数据，CUDA提供了原子操作函数——ATOM操作。

先亮出使用的例子：

[cpp] view plain copy
  1. __global__ void kernelfun()  
  2. {  
  3. __shared__ int i=0;  
  4. atomicAdd(&i, 1);  
  5. }  

如果没有加互斥机制，则同一个half warp内的线程将对i的操作混淆林乱。

用原子操作函数，可以很简单的编写自己的锁，SKD中有给出的锁结构体如下：

[cpp] view plain copy
  1. #ifndef __LOCK_H__  
  2. #define __LOCK_H__  
  3. #include "cuda_runtime.h"  
  4. #include "device_launch_parameters.h"  
  5. #include "atomic_functions.h"  
  6. struct Lock {  
  7.     int *mutex;  
  8.     Lock( void ) {  
  9.         HANDLE_ERROR( cudaMalloc( (void**)&mutex, sizeof(int) ) );  
  10.         HANDLE_ERROR( cudaMemset( mutex, 0, sizeof(int) ) );  
  11.     }  
  12.     ~Lock( void ) {  
  13.         cudaFree( mutex );  
  14.     }  
  15.     __device__ void lock( void ) {  
  16.         while( atomicCAS( mutex, 0, 1 ) != 0 );  
  17.     }  
  18.     __device__ void unlock( void ) {  
  19.         atomicExch( mutex, 0 );  
  20.     }  
  21. };  
  22. #endif  


8 CUDA软件体系结构

 

CUDA: GPU高性能运算

9 利用好现有的资源

如果连开方运算都需要自己去编写程序实现，那么我相信程序员这个职业将会缩水，没有人愿意去干这种活。我想，程序员需要学会“偷懒”，现有的资源必须学会高效率的使用。当c++出现了STL库，c++程序员的开发效率可以说倍增，而且程序稳定性更高。

CUDA有提供给我们什么了吗？给了，其实给了很多。

先介绍几个库：CUFFT、CUBLAS、CUDPP。

这里我先不详细学习这些库里到底有哪些函数，但是，大方向是需要了解的，不然找都不知道去哪儿找。CUFFT是傅里叶变换的库，CUBLAS提供了基本的矩阵和向量运算，CUDPP提供了常用的并行排序、搜索等。

CUDA4.0以上，提供了一个类似STL的模板库，初步窥探，只是一个类似vector的模板类型。有map吗？map其实是一个散列表，可以用hashtable去实现这项机制。

SDK里面有很多例子，包括一些通用的基本操作，比如InitCUDA等，都可以固化成函数组件，供新程序的调用。

具体的一些可以固化的东西，我将在以后的学习中归纳总结，丰富自己的CUDA库！

 http://blog.csdn.net/huangfengxiao/article/details/8732789

http://blog.csdn.net/huangfengxiao/article/details/8732790

http://blog.csdn.net/huangfengxiao/article/details/8732791
标签：

相关文章