Nvidia Summer Camp Day2 个人心得

这是个人第二篇心得体会，先总结一下CUDA中的一些概念

Thread：线程  CUDA中所有的线程执行相同的核函数

Thread Block：多个thread组成一个block，同一个block中threads可以使用_syncthreads()同步，也可以通过shared memory通信。执行在一个SM中。

Thread Grid：多个block组成一个grid，grid中的blocks可能分布在多个SM中执行。

kernel_function<<<dim_Grid，dim_Block，Ns，S >>>：执行配置运算符，用来传递内核函数的执行参数。执行配置有四个参数。

第一个参数dim_Grid声明一个grid中有多少block。

第二个参数dim_Block声明一个block中有多少thread。

第三个参数Ns是一个可选参数，用于设置每个block除了静态分配的shared Memory以外，最多能动态分配的shared memory大小，单位为byte。不需要动态分配时该值为0。

最后一个参数S是一个cudaStream_t类型的可选参数，初始值为零，表示该核函数处在哪个流之中。

CUDA的执行流程：加载核函数，将Grid分配到一个Device，根据<<<>>>内的执行设置的第一个参数，将Block分配到SM中。一个Block内的线程一定会在同一个SM内，一个SM可以有很多个Block。根据<<<>>>内的执行设置的第二个参数，Warp调度器会调用线程。Warp调度器为了提高运行效率，会将每32个线程分为一组，称作一个warp。

举个简单的例子

__global__ void add( int *a, int *b, int *c ) 
{
    c[threadIdx.x] = a[threadIdx.x] + b[threadIdx.x];
}
    add<<<1,3>>>( a, b, c);

数组相加的代码中<<<1,4>>>，1代表1个block，3代表3个thread，将1个block分配到SM中，Warp调度器调用3个thread，在device中运行的样子为：

thread0：c[0] = a[0] + b[0];

thread1：c[1] = a[1] + b[1];

thread2：c[2] = a[2] + b[2];

接下来说一下CUDA中的线程索引：

 

 如何计算图中5 的index呢？很显然我们只需要将blockDim.x *blockIdx.x + threadidx.x即可得到index了。

下面通过数组相加的代码来具体分析

//声明一个核函数
void __global__ add(const double *x, const double *y, double *z, int count)
{
    const int n = blockDim.x * blockIdx.x + threadIdx.x;//计算索引
	if( n < count)
	{
	    z[n] = x[n] + y[n];
	}
}

    //host部分

    const int N = 1000;
    const int M = sizeof(double) * N;
    double *h_x = (double*) malloc(M);//host中开辟内存空间
    double *h_y = (double*) malloc(M);
    double *h_z = (double*) malloc(M);

    for (int n = 0; n < N; ++n)//初始化数组
    {
        h_x[n] = 1;
        h_y[n] = 2;
    }

    //device部分

    double *d_x, *d_y, *d_z;
    cudaMalloc((void **)&d_x, M);//device中开辟内存空间
    cudaMalloc((void **)&d_y, M);
    cudaMalloc((void **)&d_z, M);
    cudaMemcpy(d_x, h_x, M, cudaMemcpyHostToDevice);//数据从Host传输到Device（HostToDevice）
    cudaMemcpy(d_y, h_y, M, cudaMemcpyHostToDevice);

    //设置block_size grid_size

    const int block_size = 128;
    const int grid_size = (N + block_size - 1) / block_size;
    add<<<grid_size, block_size>>>(d_x, d_y, d_z, N);

    cudaMemcpy(h_z, d_z, M, cudaMemcpyDeviceToHost);//计算完的结果传输回Host（DeviceToHost）


    

（关于void**：d_x指向GPU，用来保存/写入，为了在GPU中做数据处理，我们希望这个指针指向显存。GPU和CPU之间通过PCI-E相连。但是二者储存的地址不同。由于指令在CPU中处理，就必须让CPU知道在GPU中的地址位置。）
最后记得free内存空间！！！

如果我们的数据过大，thread不够用怎么办？

那就让它“一心多用”

code：

int index = blockDim.x * blockIdx.x + threadIdx.x;
	int stride = blockDim.x * gridDim.x;
	for(; index <n; index +=stride)
		z[index] = x[index] + y[index];

重点！接下来开始矩阵乘法的GPU实现：

首先我们先看CPU中的矩阵乘法该如何实现，code：

void cpu_matrix_mult(int *h_a, int *h_b, int *h_result, int m, int n, int k) 
{
    for (int i = 0; i < m; ++i) 
    {
        for (int j = 0; j < k; ++j) 
        {
            int tmp = 0.0;
            for (int h = 0; h < n; ++h) 
            {
                tmp += h_a[i * n + h] * h_b[h * k + j];
            }
            h_result[i * k + j] = tmp;
        }
    }
}

大佬们应该对此代码一目了然，由于本人代码写的较少，对于此部分还是花了一些时间来理解，我将我的个人理解写出来：

我们先定义两个矩阵 h_a（m * n）与 h_b（n * k），h_result（m * k）来保存结果。

计算两个矩阵的索引：首先是 a 矩阵的，m 行n 列，i 表示所在的行数，n 为定值，i * n 表示 h 之前的所有行有多少索引，i * n + h 就变得很好理解了，后面的 h* k + j 同理。tmp 用来储存 result。

可以看出，在CPU中，用for循环会一个接着一个处理数值，直到处理完所有的线程。下面我们来看在GPU中我们如何优化矩阵乘法。

blockIdx（1，1）所以 blockIdx.x = 1，blockIdx.y = 1 ；

blockDim （4，3）所以blockDim.x = 4，blockDim.y = 3 ；

threadIdx（2，0）所以threadIdx.x = 2，threadIdx.y = 0；

得到Thread_x = 1 * 4 + 2， Thread_y = 1 * 0 + 3 = 3。

code：

__global__ void gpu_matrix_mult(int *a,int *b, int *c, int m, int n, int k)
{ 
    int row = blockIdx.y * blockDim.y + threadIdx.y; //行
    int col = blockIdx.x * blockDim.x + threadIdx.x; //列
    int sum = 0;
    if( col < k && row < m) 
    {
        for(int i = 0; i < n; i++) 
        {
            sum += a[row * n + i] * b[i * k + col];  //与CPU原理相同
        }
        c[row * k + col] = sum;
    }
} 

到这里我们已经完成了kernel。

函数的实现：

    int m=;
    int n=;
    int k=;//自行赋值

    int *h_a, *h_b, *h_c, *h_cc;//（h代表host）
    cudaMallocHost((void **) &h_a, sizeof(int)*m*n);
    cudaMallocHost((void **) &h_b, sizeof(int)*n*k);
    cudaMallocHost((void **) &h_c, sizeof(int)*m*k);
    cudaMallocHost((void **) &h_cc, sizeof(int)*m*k);//CPU开辟内存空间


    //随机数初始化矩阵
    for (int i = 0; i < m; ++i) {
        for (int j = 0; j < n; ++j) {
            h_a[i * n + j] = rand() % 1024;
        }
    }

    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < k; ++j) {
            h_b[i * k + j] = rand() % 1024;
        }
    }

    int *d_a, *d_b, *d_c;//（d代表device）
    cudaMalloc((void **) &d_a, sizeof(int)*m*n);
    cudaMalloc((void **) &d_b, sizeof(int)*n*k);
    cudaMalloc((void **) &d_c, sizeof(int)*m*k);//GPU开辟内存空间

    // 将矩阵A与矩阵B从Host传到Device
    cudaMemcpy(d_a, h_a, sizeof(int)*m*n, cudaMemcpyHostToDevice);
    cudaMemcpy(d_b, h_b, sizeof(int)*n*k, cudaMemcpyHostToDevice);

    unsigned int grid_rows = (m + BLOCK_SIZE - 1) / BLOCK_SIZE;
    unsigned int grid_cols = (k + BLOCK_SIZE - 1) / BLOCK_SIZE;
    dim3 dimGrid(grid_cols, grid_rows);
    dim3 dimBlock(BLOCK_SIZE, BLOCK_SIZE);
   
    //！！！！！
    gpu_matrix_mult<<<dimGrid, dimBlock>>>(d_a, d_b, d_c, m, n, k);    
    //将计算结果从Device传回Host
    cudaMemcpy(h_c, d_c, sizeof(int)*m*k, cudaMemcpyDeviceToHost);

最后别忘了 free ! 

可以看出，CUDA的kernel中，我们不需要再用for循环一个接着一个处理数值，类似于数组相加，每个线程，读取对应的行列数据，相乘求和。GPU并行计算相比较CPU，它可以分配多个线程，每个线程并行，所以可以达到加速的目的。

第二篇心得好像比第一篇有所进步，下一篇会总结有关 shared memory 以及原子操作的内容，欢迎大家指出问题以及交流！谢谢观看！