3、CUDA中的向量运算
3.1、两个向量加法程序
GPU main函数具有本章第一节所述的已知结构：
1、先是定义CPU和GPU上的数组和指针。设备指针指向通过cudaMalloc分配的显存。
2、然后通过cudaMemcpy函数，将前两个数组，从主机内存传输到设备显存。
3、内核启动的时候，将这些设备指针作为参数传递给它。你看到的内核启动符号（<<<>>>）里面的N和1，它们分别是启动N个块，每个块里面只有1个线程。
4、再然后通过cudaMemcpy，将内核的计算结果从设备显存传输到主机内存，注意最后这次传输的方向是反的，设备到主机。
5、最后，用cudaFree释放掉3段在显存上分配的缓冲区。
直接上代码：
#include "stdio.h" #include<iostream> #include <cuda.h> #include <cuda_runtime.h> //定义数组元素个数 #define N5 //定义GPU设备的Kernel函数 __global__ void gpuAdd(int *d_a, int *d_b, int *d_c) { int tid = blockIdx.x;// 获取数据的index if (tid < N) d_c[tid] = d_a[tid] + d_b[tid]; } int main(void) { //定义主机数组 int h_a[N], h_b[N], h_c[N]; //定义GPU设备指针 int *d_a, *d_b, *d_c; // 申请内存 cudaMalloc((void**)&d_a, N * sizeof(int)); cudaMalloc((void**)&d_b, N * sizeof(int)); cudaMalloc((void**)&d_c, N * sizeof(int)); //初始化数组 for (int i = 0; i < N; i++) { h_a[i] = 2 * i*i; h_b[i] = i; } // 把CPU下的数组复制到GPU设备 cudaMemcpy(d_a, h_a, N * sizeof(int), cudaMemcpyHostToDevice); cudaMemcpy(d_b, h_b, N * sizeof(int), cudaMemcpyHostToDevice); //确定要使用的块和线程数目 gpuAdd << <N, 1 >> > (d_a, d_b, d_c); //把GPU运输的结果复制回CPU主机上进行打印处理 cudaMemcpy(h_c, d_c, N * sizeof(int), cudaMemcpyDeviceToHost); printf("Vector addition on GPU \n"); //在控制台打印结果 for (int i = 0; i < N; i++) { printf("The sum of %d element is %d + %d = %d\n", i, h_a[i], h_b[i], h_c[i]); } //最后别忘了释放内存 cudaFree(d_a); cudaFree(d_b); cudaFree(d_c); return 0; }
所有CUDA程序都遵循与前面相同的模式。我们并行启动N个块。这意味着我们同时启动了N个执行该内核代码的线程副本。你可以通过一个现实的例子来理解这一点：假设你想把5个大盒子从一个地方转移到另一个地方。在第一种方法中，你可以通过雇佣一个人将一个盒子从一个地方带到另一个地方，然后重复5次来完成这项任务。这个方式需要时间，它类似于如何在CPU上做向量相加。现在，假设你雇了5个人，每个人都带着一个盒子。他们每个人也知道他们所携带的箱子的ID。这个方法将比前一个快得多。每个人只需要被告知他们必须携带一个带有特定ID的盒子从一个地方到另一个地方。
从上述例子里我们知道内核代码的编写方式，以及这些线程是如何在GPU上执行的。每个线程可以通过blockIdx.x内置变量来知道自己的ID。然后每个线程通过ID来索引数组，计算每对元素加法。这样，多个线程的并行计算，明显减少了数组整体的处理时间。所以，用这种并行的计算方式，比CPU上的串行计算，提高了吞吐率。
这个程序还含有额外的一个CUDA函数调用：cudaDeviceSynchronize()。为何要加这句？这是因为启动内核是一个异步操作，只要发布了内核启动命令，不等内核执行完成，控制权就会立刻返回给调用内核的CPU线程。在上述的代码中，CPU线程返回，继续执行的下一句是printf()。而再之后，在内核完成之前，进程就会结束，终止控制台窗口。所以，如果不加上这句同步函数，你就看不到任何的内核执行结果输出。在程序退出后内核生成的输出结果，将没有地方可去，你没法看到它们，因此，如果我们不包含这个指令，你将不会看到任何内核执行的printf语句的输出结果。要能看到内核生成的输出结果，我们必须包含这句同步函数。这样，内核的结果将通过可用的标准输出显示，而应用程序则会在内核执行完成之后才退出。
3.2、对比CPU代码和GPU代码的延迟
#include "stdio.h" #include<iostream> #include <cuda.h> #include <cuda_runtime.h> //定义数组元素个数 #define N10000000 //定义GPU设备的Kernel函数 __global__ void gpuAdd(int *d_a, int *d_b, int *d_c) { int tid = blockIdx.x;// 获取数据的index if (tid < N) d_c[tid] = d_a[tid] + d_b[tid]; } //定义CPU设备的函数 void cpuAdd(int *h_a, int *h_b, int *h_c) { int tid = 0; while (tid < N) { h_c[tid] = h_a[tid] + h_b[tid]; tid += 1; } } int main(void) { //定义主机数组 int h_a[N], h_b[N], h_c[N]; //定义GPU设备指针 int *d_a, *d_b, *d_c; // 申请内存 cudaMalloc((void**)&d_a, N * sizeof(int)); cudaMalloc((void**)&d_b, N * sizeof(int)); cudaMalloc((void**)&d_c, N * sizeof(int)); //初始化数组 for (int i = 0; i < N; i++) { h_a[i] = 2 * i*i; h_b[i] = i; } /*****************************************CPU运算****************************************/ clock_t start_h = clock(); //Calling CPU function for vector addition cpuAdd(h_a, h_b, h_c); clock_t end_h = clock(); double time_h = (double)(end_h - start_h) / CLOCKS_PER_SEC; //printf("Vector addition on CPU\n"); //for (int i = 0; i < N; i++) { //printf("The sum of %d element is %d + %d = %d\n", i, h_a[i], h_b[i], h_c[i]); //} /*****************************************GPU运算****************************************/ // 把CPU下的数组复制到GPU设备 cudaMemcpy(d_a, h_a, N * sizeof(int), cudaMemcpyHostToDevice); cudaMemcpy(d_b, h_b, N * sizeof(int), cudaMemcpyHostToDevice); clock_t start_d = clock(); //确定要使用的块和线程数目并在GPU处理运行 gpuAdd << <N, 1 >> > (d_a, d_b, d_c); cudaThreadSynchronize(); clock_t end_d = clock(); double time_d = (double)(end_d - start_d) / CLOCKS_PER_SEC; printf("No of Elements in Array:%d\n Device time %f Seconds \n host time %f Seconds\n", N, time_d, time_h); //把GPU运输的结果复制回CPU主机上进行打印处理 cudaMemcpy(h_c, d_c, N * sizeof(int), cudaMemcpyDeviceToHost); //printf("Vector addition on GPU \n"); //在控制台打印结果 //for (int i = 0; i < N; i++) { //printf("The sum of %d element is %d + %d = %d\n", i, h_a[i], h_b[i], h_c[i]); //} //最后别忘了释放显存 cudaFree(d_a); cudaFree(d_b); cudaFree(d_c); return 0; }
我们可以从输出结果看到，执行时间或叫吞吐量，从25ms提高到在GPU上几乎1ms实现相同的功能。这证明了我们之前在理论上提到的GPU并行执行代码有助于提高吞吐量。
3.3、并行通信模式
当多个线程并行执行时，它们遵循一定的通信模式，指导它们在显存里哪里输入，哪里输出。我们将讨论每个通信模式。它将帮助你识别通信模式相关的应用程序，以及如何编写代码。
1、映射模式（一对一操作：基本运算使用）
在这种通信模式中，每个线程或任务读取单一输入，产生一个输出。基本上，它是一个一对一的操作。前面部分中向量加法程序和向量元素平方的程序，就是Map模式的例子。Map的代码模式看起来如下：
2、收集（多对一操作，移动平均使用）
在此模式中，每个线程或者任务，具有多个输入，并产生单一输出，保存到存储器的单一位置中。假设你想写一个求3数据的MA操作的程序。这就是一个Gather操作的例子。每个线程读取3个输入数据，产生单一的结果数据保存到显存。因此，在输入端有数据复用。它基本上是一个多对一的操作。Gather模式的代码看起来如下：
3、分散式（一对多操作，按序排列使用）
在Scatter模式中，线程或者任务读取单一输入，但向存储器产生多个输出。数组排序就是一个Scatter操作的例子。它也可以叫作1对多操作。Scatter模式的代码看起来如下：
4、蒙板（图像使用）
当线程或者任务要从数组中读取固定形状的相邻元素时，这叫stencil模式。这种模式在图像处理中非常有用，例如当你想用一个3×3或者5×5的（滤波）窗口时（对整个图像进行滑动处理的时候）。它是Gather操作的一种模式，所以代码的语法和Gather很相似。
5、转置（一对一操作，基本使用）
当原始输入矩阵是行主序的时候，如果需要输出得到一个列主序的矩阵，则应当进行转置操作。如果你有一个结构数组（SoA），而你想把它转换成一个数组结构（AoS），它是特别有用的。这也是一个一对一的操作。Transpose模式的代码看起来如下：
