CUDA编程，从零开始（2）
本次记录我自己写的矩阵乘法的cuda程序，以及CPU，GPU加速以及GPU封装库函数的运行时间对比。
这里面还涉及到一些小问题，暂时写到注释中，以后有时间会再整理出一篇来单独说。
记录我的笔记本版本信息
Win10+Vs2015+Cuda9.1
记录我的测试矩阵信息
M1(1024,1024) * M2(1024,1024)
CPU
void matMul_C(float *mul_t, float *mul1, float *mul2)
{
for (int i = 0; i < M; i++)
{
for (int j = 0; j < M; j++)
{
float sum0 = 0;
for (int k = 0; k < M; k++)
{
sum0 = sum0 + mul1[i*M + k] * mul2[k*M + j];
}
mul_t[i*M + j] = sum0;
}
}
}
主函数部分计时：
printf("Start:\n");
start = clock();
matMul_C(mul_t, mul1, mul2);
stop = clock();
time_c = (double)(stop - start) / CLOCKS_PER_SEC;
printf("The time for CPU: %lfs\n", time_c);
GPU并行1
__global__ void matMul_G(float *mul_t, float *mul1, float *mul2)
{
int i = threadIdx.x + blockIdx.x * blockDim.x;
int j = threadIdx.y + blockIdx.y * blockDim.y;
int k;
for (k = 0; k < M; k++)
{
mul_t[i*M + j] = mul1[i*M + k] * mul2[k*M + j] + mul_t[i*M + j];
}
}
主函数部分计时与调用：
cudaEvent_t start_g, stop_g;
cudaEvent_t start_a, stop_a;
float costtime;
float costtime_g;
dim3 DimGrid(32, 32);
dim3 DimBlock(32, 32);
float *dev_t, *dev1, *dev2, *dev_G = (float*)malloc(sizeof(float)*M*M);
cudaEventCreate(&start_g);
cudaEventCreate(&stop_g);
cudaEventCreate(&start_a);
cudaEventCreate(&stop_a);
/*为GPU设备分配内存，这里的指针指向设备端地址*/
cudaMalloc((void **)&dev_t, sizeof(float)*M*M);
cudaMalloc((void **)&dev1, sizeof(float)*M*M);
cudaMalloc((void **)&dev2, sizeof(float)*M*M);
/*创建事件*/
cudaEventRecord(start_g, 0);
/*数据拷贝*/
cudaMemcpy(dev1, mul1, sizeof(float)*M*M,cudaMemcpyHostToDevice);
cudaMemcpy(dev2, mul2, sizeof(float)*M*M, cudaMemcpyHostToDevice);
cudaMemcpy(dev_t, mul_t2, sizeof(float)*M*M, cudaMemcpyHostToDevice);
/*开始计时，调用核函数*/
cudaEventRecord(start_a, 0);
matMul_G << <DimGrid, DimBlock >> > (dev_t, dev1, dev2);
cudaDeviceSynchronize();
cudaEventRecord(stop_a, 0);
cudaMemcpy(dev_G, dev_t, sizeof(float)*M*M, cudaMemcpyDeviceToHost);
cudaEventRecord(stop_g, 0);
cudaEventSynchronize(stop_g);
/*需要注意的是函数cudaEventSynchronize() 不可或缺，因为CUDA的kernel函数是以异步方式执行的，调用后立刻返回，
这会导致计时不准确。cudaEventSynchronize(stop)会使得直到GPU执行完cudaEventRecord(stop, 0)之前的所有语句时，
事件stop才会被记录下来，即起到同步的作用。*/
cudaEventElapsedTime(&costtime, start_a, stop_a);
cudaEventDestroy(start_a);
cudaEventDestroy(stop_a);
cudaEventElapsedTime(&costtime_g, start_g, stop_g);
cudaEventDestroy(start_g);
cudaEventDestroy(stop_g);
printf("The time for gpu_mul is: %lf ms\n", costtime);
printf("The time for gpu is: %lf ms\n", costtime_g);
GPU并行2——shared memory的使用
__global__ void matMul_Gs(float *mul_t, float *mul1, float *mul2)
{
/*Thread 所在 block 的 location*/
int block_x = blockIdx.x;
int block_y = blockIdx.y;
/*Thread 的location*/
int thread_x = threadIdx.x;
int thread_y = threadIdx.y;
int height_A = M, width_A = M, width_B = M;
/*A矩阵子矩阵的起始下标*/
int begin_a = block_y * blockDim.y * width_A;
/*A矩阵子矩阵的终止下标（就是A矩阵一次运算一行，对应着B矩阵一次运算一列）*/
int end_a = begin_a + width_A - 1;
/*A矩阵下标一次移动的步长，一次处理一个子矩阵*/
int step_a = blockDim.x;
/*B 矩阵子矩阵对应的起始下标*/
int begin_b = block_x * blockDim.x;
/*B 矩阵子矩阵对应的步长*/
int step_b = blockDim.y * width_B;
int result_temp = 0;
/*循环次数等于widthA / bs，把长向量点积运算转化为短向量点积后的和*/
for (int index_a = begin_a, int index_b = begin_b; index_a < end_a; index_a += step_a, index_b += step_b)
{
__shared__ int A[BLOCK_SIZE][BLOCK_SIZE];
__shared__ int B[BLOCK_SIZE][BLOCK_SIZE];
A[thread_y][thread_x] = mul1[index_a + thread_y * width_A + thread_x];
B[thread_y][thread_x] = mul2[index_b + thread_y * width_B + thread_x];
__syncthreads();
for (int i = 0; i < BLOCK_SIZE; i++)
{
result_temp += A[thread_y][i] * B[i][thread_x];
}
__syncthreads();
}
int begin_result = block_y * blockDim.y * width_B + begin_b;
mul_t[begin_result + thread_y * width_B + thread_x] = result_temp;
}
主函数部分计时与调用：
cudaEvent_t start_g2, stop_g2;
cudaEvent_t start_a2, stop_a2;
float costtime2;
float costtime_g2;
dim3 DimGrid2(32, 32);
dim3 DimBlock2(32, 32);
float *dev_t2, *dev12, *dev22, *dev_G2 = (float*)malloc(sizeof(float)*M*M);
cudaEventCreate(&start_g2);
cudaEventCreate(&stop_g2);
cudaEventCreate(&start_a2);
cudaEventCreate(&stop_a2);
/*为GPU设备分配内存，这里的指针指向设备端地址*/
cudaMalloc((void **)&dev_t2, sizeof(float)*M*M);
cudaMalloc((void **)&dev12, sizeof(float)*M*M);
cudaMalloc((void **)&dev22, sizeof(float)*M*M);
/*为什么第一个参数是二级指针呢？
因为devPtrA是一个CPU上的指针，我现在分配了一块GPU内存空间，假设这个空间的首地址是 #0024，
我想要修改devPtrA的值怎么办呢？C语言的基础告诉我们，函数想要修改一个变量的值，
就得把他的指针传进去，而不是仅仅作为形参传过去。所以这里是二级指针，
因为我要把一级指针的值修改为 #0024。
两个函数的最大的不同是函数的返回值类型，malloc返回的是void *，而cudaMalloc返回的是cudaError_t
（它不是指针类型哦，对其实根本原因就是这个），CUDA中的函数好像都是采用这种结构返回值类型为cudaError_t，
但现在我们是想得到函数cudaMalloc在device上申请的内存，这可咋搞？
在函数中为形参赋值是不会在实参中繁盛变化的，但是指针传递的是地址，我们操作了某个地址的数据，
实际上是真的改变了指定地址的数据。像这个申请显存的函数，第一个参数传递的是devPtrA这个指针的地址，
然后改变这个地址的内容就会带给实参真正的改变。既然要把申请的内存地址返回给一个指针，
那就用二重指针解决吧，这样既解决了CUDA函数返回类型的统一，也解决了申请的内存地址返回给一个指针问题
*/
/*创建事件*/
cudaEventRecord(start_g2, 0);
/*数据拷贝*/
cudaMemcpy(dev12, mul1, sizeof(float)*M*M, cudaMemcpyHostToDevice);
cudaMemcpy(dev22, mul2, sizeof(float)*M*M, cudaMemcpyHostToDevice);
cudaMemcpy(dev_t2, mul_t2, sizeof(float)*M*M, cudaMemcpyHostToDevice);
/*开始计时，调用核函数*/
cudaEventRecord(start_a2, 0);
matMul_Gs << <DimGrid2, DimBlock2 >> > (dev_t2, dev12, dev22);
cudaDeviceSynchronize();
cudaEventRecord(stop_a2, 0);
cudaMemcpy(dev_G2, dev_t2, sizeof(float)*M*M, cudaMemcpyDeviceToHost);
cudaEventRecord(stop_g2, 0);
cudaEventSynchronize(stop_g2);
/*需要注意的是函数cudaEventSynchronize() 不可或缺，因为CUDA的kernel函数是以异步方式执行的，调用后立刻返回，
这会导致计时不准确。cudaEventSynchronize(stop)会使得直到GPU执行完cudaEventRecord(stop, 0)之前的所有语句时，
事件stop才会被记录下来，即起到同步的作用。*/
cudaEventElapsedTime(&costtime2, start_a2, stop_a2);
cudaEventDestroy(start_a2);
cudaEventDestroy(stop_a2);
cudaEventElapsedTime(&costtime_g2, start_g2, stop_g2);
cudaEventDestroy(start_g2);
cudaEventDestroy(stop_g2);
printf("The time for gpu_mul_S is: %lf ms\n", costtime2);
printf("The time for gpu_S is: %lf ms\n", costtime_g2);
GPU并行3——cublas库函数调用
- 添加头文件 #include <cublas_v2.h>
- 右击项目–>点击属性–>选择链接器–>输入–>附加依赖项–>cublas.lib
主函数部分调用：
int N = 1024;
cudaEvent_t start_g3, stop_g3;
cudaEvent_t start_a3, stop_a3;
float costtime3;
float costtime_g3;
float *dev_t3, *dev13, *dev23, *dev_G3 = (float*)malloc(sizeof(float)*M*M);
float alpha = 1.0;
float beta = 0.0;//C=α∗A∗B+β∗C
cudaEventCreate(&start_g3);
cudaEventCreate(&stop_g3);
cudaEventCreate(&start_a3);
cudaEventCreate(&stop_a3);
/*为GPU设备分配内存，这里的指针指向设备端地址*/
cudaMalloc((void **)&dev_t3, sizeof(float)*M*M);
cudaMalloc((void **)&dev13, sizeof(float)*M*M);
cudaMalloc((void **)&dev23, sizeof(float)*M*M);
cublasHandle_t handle;//管理一个cublas上下文的句柄
cublasCreate_v2 (&handle);//初始化
/*创建事件*/
cudaEventRecord(start_g3, 0);
/*数据拷贝*/
cudaMemcpy(dev13, mul1, sizeof(float)*M*M, cudaMemcpyHostToDevice);
cudaMemcpy(dev23, mul2, sizeof(float)*M*M, cudaMemcpyHostToDevice);
cudaMemcpy(dev_t3, mul_t2, sizeof(float)*M*M, cudaMemcpyHostToDevice);
//cudaMemset
cudaEventRecord(start_a3, 0);
cublasSgemm_v2 (handle, CUBLAS_OP_N, CUBLAS_OP_N, N, N, N,
&alpha, dev23, N, dev13, N, &beta, dev_t3, N);//C=A*B
/*可以等于CUBLAS_OP_T表示需要转置后运算，或者不需要CUBLAS_OP_N；
m : mul1有多少行； //矩阵B的列数
n : mul2有多少列； //矩阵A的行数
k : mul1的列数，也是mul2的行数；//矩阵A的列数
lda是A的leading dimension ，既然是列优先，那么就是A的行数（A的一列有多少个元素）；ldb同理*/
cudaDeviceSynchronize();
cudaEventRecord(stop_a3, 0);
cudaMemcpy(dev_G3, dev_t3, sizeof(float)*M*M, cudaMemcpyDeviceToHost);
cudaDeviceSynchronize();
cudaEventRecord(stop_g3, 0);
cudaEventSynchronize(stop_g3);
/*需要注意的是函数cudaEventSynchronize() 不可或缺，因为CUDA的kernel函数是以异步方式执行的，调用后立刻返回，
这会导致计时不准确。cudaEventSynchronize(stop)会使得直到GPU执行完cudaEventRecord(stop, 0)之前的所有语句时，
事件stop才会被记录下来，即起到同步的作用。*/
cublasDestroy_v2 (handle);
cudaEventElapsedTime(&costtime3, start_a3, stop_a3);
cudaEventDestroy(start_a3);
cudaEventDestroy(stop_a3);
cudaEventElapsedTime(&costtime_g3, start_g3, stop_g3);
cudaEventDestroy(start_g3);
cudaEventDestroy(stop_g3);
printf("The time for gpu_mul cublas is: %lf ms\n", costtime3);
printf("The time for gpu cublas is: %lf ms\n", costtime_g3);
最后一定记得释放内存
cudaFree(dev1);
cudaFree(dev2);
cudaFree(dev_t2);
cudaFree(dev12);
cudaFree(dev22);
cudaFree(dev_t3);
cudaFree(dev13);
cudaFree(dev23);
free(mul_t);
free(mul1);
free(mul2);
free(mul_t2);
free(dev_G);
free(dev_G2);
free(dev_G3);
计时对比结果
1278
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
