一、使用全局内存
矩阵乘法,即用矩阵A每行与矩阵B的每列,依次做乘积累加就可以得到各个元素的值。在CPU上用三层循环实现。这里是将二维数组用一维的形式表示,即按行存储。
size_t size = WIDTH*WIDTH * sizeof(int); int *h_A = (int *)malloc(size); int *h_B = (int*)malloc(size); int *h_C = (int*)malloc(size); int *h_d_Answer = (int *)malloc(size); if (h_A == NULL || h_B == NULL || h_C == NULL || h_d_Answer == NULL) { exit(EXIT_FAILURE); } for (int i = 0; i < WIDTH; i++) { for (int j = 0; j < WIDTH; j++) { h_A[i*WIDTH + j] = i*WIDTH + j; h_B[i*WIDTH + j] = 1; } } for (int i = 0; i < WIDTH; i++) { for (int j = 0; j < WIDTH; j++) { int sum = 0; for (int k = 0; k < WIDTH; k++) { sum += h_A[i*WIDTH + k] * h_B[k*WIDTH + j]; } h_C[i*WIDTH + j] = sum; } } printf("CPU answer:\n"); 使用全局内存在GPU上的作法是,用idx idy作为二维数组的索引形式,这样就可以通过[idy*WIDTh + idx]确定到具体的元素,每个线程读取A的一行和B的一列,然后一个循环作累加。在GPU对应的kernel 函数:
__global__ void kernelMatrix(int *A, int *B, int *C){ int idx = threadIdx.x + blockDim.x*blockIdx.x;//col number int idy = threadIdx.y + blockDim.y*blockIdx.y;//row number if (idx < WIDTH && idy < WIDTH) { int sum = 0; for (int k = 0; k < WIDTH; k++) { sum += A[idy*WIDTH + k] * B[k*WIDTH + idx]; } C[idy*WIDTH + idx] = sum; } } 二、使用共享内存每个Block计算一个方阵的子矩阵,大小为BLOCKDIM,由BLOCK的共享内存装载数据。
线程配置示意图:
数据装载过程示意图:
======================================================================
• 外部 for 循环控制所有的线程读取的次数WIDTH/TILE_WIDTH • 内部 for 循环进行共享内存数据的相乘并累加实现代码如下:
__global__ void kernelMatrixShare(int *A, int *B, int *C) { __shared__ int ds_M[BLOCKDIM][BLOCKDIM]; __shared__ int ds_N[BLOCKDIM][BLOCKDIM]; int idx = threadIdx.x + BLOCKDIM*blockIdx.x; int idy = threadIdx.y + BLOCKDIM*blockIdx.y; int tx = threadIdx.x; int ty = threadIdx.y; int sum = 0; for (int m = 0; m < WIDTH / BLOCKDIM; m++) { ds_M[ty][tx] = A[idy*WIDTH + m*BLOCKDIM+tx];//A[idy][m*BLOCKDIM + tx]; ds_N[ty][tx] = B[idx + (m*BLOCKDIM + ty)* WIDTH]; __syncthreads(); for (int k = 0; k < BLOCKDIM; k++) { sum += ds_M[ty][k] * ds_N[k][tx]; } __syncthreads(); } //get one value C[idy*WIDTH + idx] = sum; }实验对比结果: