CUDA中的向量计算与并行通信模式

• 本节开始，我们将利用GPU的并行能力，对其执行向量和数组操作
• 讨论每个通信模式，将帮助你识别通信模式相关的应用程序，以及如何编写代码

1.两个向量加法程序

• 先写一个通过cpu实现向量加法的程序
• 如下所示，向量相加实际上是模仿GPU的写法，在GPU中，tid 代表特定的某个线程的ID。
• 如果你的cpu是双核的，可以在每个核心上运行一个线程，分别将tid初始化为0和1，然后每次循环的时候+2，这样的话可以实现一个核激素那偶数元素的和，一个核计算基数元素的和，通过两个线程的实现并行计算

#include "stdio.h"
#include<iostream>
//Defining Number of elements in Array
#define N	5
//Defining vector addition function for 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];//Initializing two arrays for additionfor (int i = 0; i < N; i++) {h_a[i] = 2 * i*i;h_b[i] = i;}//Calling CPU function for vector additioncpuAdd (h_a, h_b, h_c);//Printing Answerprintf("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]);}return 0;
}

• 而总所周知，NVIDIA GPU包含多个块，每个块又包含多个线程，因此可以通过GPU实现更多线程并行计算向量的和，最大程度提高速度
• 可以将代码修改为核函数如下：

#include "stdio.h"
#include<iostream>
#include <cuda.h>
#include <cuda_runtime.h>//Defining number of elements in Array
#define N	5//Defining Kernel function for vector addition
__global__ void gpuAdd(int* d_a, int* d_b, int* d_c) {//Getting block index of current kernelint tid = blockIdx.x;	// handle the data at this indexif (tid < N)d_c[tid] = d_a[tid] + d_b[tid];
}int main(void) {//定义主机数组变量int h_a[N], h_b[N], h_c[N];//定义设备指针变量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));//Initializing Arraysfor (int i = 0; i < N; i++) {h_a[i] = 2 * i * i;h_b[i] = i;}// Copy input arrays from host to device memorycudaMemcpy(d_a, h_a, N * sizeof(int), cudaMemcpyHostToDevice);cudaMemcpy(d_b, h_b, N * sizeof(int), cudaMemcpyHostToDevice);//设置内核参数为5个块，每个块一个线程 ，并像核函数传递参数gpuAdd << <N, 1 >> > (d_a, d_b, d_c);//将计算结果从显卡拷贝到主机cudaMemcpy(h_c, d_c, N * sizeof(int), cudaMemcpyDeviceToHost);printf("Vector addition on GPU \n");//Printing result on consolefor (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]);}//Free up memorycudaFree(d_a);cudaFree(d_b);cudaFree(d_c);return 0;
}

• 以上可以发现，通过GPU并行运算，或者多个线程的并行计算，明显的减少了数组的处理时间。比起CPU上的串行计算，提高了吞吐率
• 在此说一下吞吐量的含义：只对网络、设备端口、虚电路或者其他设施，单位时间内成功的传送数据的数量（以比特、字节、分贝等测量）

2. 对比CPU代码和GPU代码的延迟

• CPU的加法程序和GPU的加法程序都是以一个模块化的方式来编写的
• N的值较小时，看不出cpu与GPU的差异，但是当N值很大时，会发现两者计算效率的显著差异
• 下边将展示如何为并行计算计时并对两者时间进行比较

clock_t start_d = clock();printf("Doing GPU Vector add\n");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 second \n Host time %f second \n ",N, time_d, time_h);

3. 对向量的每个元素进行平方

• 前边调用内核函数时启用了N个块，每个块一个线程执行计算；另一种也可以只启动1个块，块里边有N个线程，淡然也可以启用N个块，每个块M个线程
• 下边通过启用一个块中的N个线程来执行向量每个元素的平方运算

#include "stdio.h"
#include<iostream>
#include <cuda.h>
#include <cuda_runtime.h>
//Defining number of elements in Array
#define N	5
//Kernel function for squaring number
__global__ void gpuSquare(float *d_in, float *d_out) {//Getting thread index for current kernelint tid = threadIdx.x;	// handle the data at this indexfloat temp = d_in[tid];d_out[tid] = temp*temp;
}int main(void) {//Defining Arrays for hostfloat h_in[N], h_out[N];//Defining Pointers for devicefloat *d_in, *d_out;// allocate the memory on the gpucudaMalloc((void**)&d_in, N * sizeof(float));cudaMalloc((void**)&d_out, N * sizeof(float));//Initializing Arrayfor (int i = 0; i < N; i++) {h_in[i] = i;}//Copy Array from host to devicecudaMemcpy(d_in, h_in, N * sizeof(float), cudaMemcpyHostToDevice);//Calling square kernel with one block and N threads per blockgpuSquare << <1, N >> >(d_in, d_out);//Coping result back to host from device memorycudaMemcpy(h_out, d_out, N * sizeof(float), cudaMemcpyDeviceToHost);//Printing result on consoleprintf("Square of Number on GPU \n");for (int i = 0; i < N; i++) {printf("The square of %f is %f\n", h_in[i], h_out[i]);}//Free up memorycudaFree(d_in);cudaFree(d_out);return 0;
}

• 需注意

  • 每当使用这种方式启动N个线程并行的时候，需要注意每个块的最大线程不超过 512 或 1024
  • 现在所有计算能力/显卡算力在 3.0 - 7.5 的GPU卡，每个块最大1024个线程
  • 如果N是2000，而你的GPU卡线程的最大数量是512，那么不能写成 << <12 000 > >>，而应该使用<< <4,500 > >>，应该理性的选择合适数量的块和每个块具有的线程数量

4. 并行通信模式

• 当多个线程并行执行时，它们遵循一定的通信模式，知道它们在显存里哪里输入，哪里输出

4.1 映射

• 一对一操作，每个线程或任务读取单一输入，产生单一输出，就是Map模式

d_out[i] = d_in[i] * 2

4.2 收集

• 此模式下，每个线程或者任务，具有多个输入，并产生单个输出，保存到存储器的单一位置，即Gather模式：

out[i] = (in[i-1] + in[i] + in[i+1]) / 3

4.3 分散式

• Scatter 模式，线程或者任务读取单一输入，单项存储器产生多个输出，比如数组排序：

out[i-1] += 2 * in[i] and out[i+1] += 3 * in[i]

4.4 蒙版

• 当线程或者任务要从数组中读取固定形状的相邻元素时，这叫stencil模式，在图像处理中非常有用。比如想用一个3X3或者5X5的窗口进行滑动滤波
• 代码类似Gather

4.5 转置

• 当想要输入矩阵行主序，输出矩阵想要列主序，或者有一个结构数组（SoA)，想转换成一个数组结构(AoS)，它是特别有用的。Transpose模式如下：

out[i+j*128] = in[j + i*128]