CUDA编程实践-3

矩阵转置求解

实现了两个版本的转置。一种是没有进行内存优化，直接上gpu的版本。

一种是考虑访存优化，利用share-memory进行优化。

先把host端的代码粘出来：

host端代码：

#include<stdio.h>
#include <stdlib.h>

#include <cuda_runtime.h>
#include <helper_cuda.h>
#define BLOCK 16

void checkMalloc(float *ptr){
if(ptr ==  NULL){
printf("malloc failed\n");
exit(0);
}
}

void dealError(cudaError_t err){
if(err != cudaSuccess){
printf("%s\n",cudaGetErrorString(cudaGetLastError()));
}

}

void initFloatVec(float *h_src, int row, int col){
int i, j;
for(i=0; i<row; i++)
for(j=0; j<col; j++){
h_src[i + j * row] = i;
}

}
void printFloatVec(float *h_src, int row, int col){
int i, j;
for(i=0; i<row; i++){
for(j=0; j<col; j++)
printf("%f ", h_src[i + j * row]);
printf("\n");
}
printf("\n\n");
}

int main(){
//int row = 10, col = 8;
int row = 4, col = 2;
//  int row = 8, col = 8;
int memSize = row*col*sizeof(float);

cudaError_t err = cudaSuccess;

float *h_src =  NULL;
float *h_dst =  NULL;

float *d_src =  NULL;
float *d_dst =  NULL;

dim3 threadsPerBlock(BLOCK, BLOCK);
dim3 blocksPerGrid((row+BLOCK-1)/BLOCK, (col+BLOCK-1)/BLOCK);
/* malloc */
h_src = (float*)malloc(memSize);
h_dst = (float*)malloc(memSize);
checkMalloc(h_src);
checkMalloc(h_dst);

err = cudaMalloc((void**)&d_src, memSize);
dealError(err);

err = cudaMalloc((void**)&d_dst, memSize);
dealError(err);
/*init h_src*/
initFloatVec(h_src, row, col);
printFloatVec(h_src, row, col);
/*transfer to device memory*/
err = cudaMemcpy(d_src, h_src, memSize, cudaMemcpyHostToDevice);
dealError(err);

transposeOpt<<<blocksPerGrid, threadsPerBlock>>>(d_src, d_dst, row, col);

err = cudaMemcpy(h_dst, d_dst, memSize, cudaMemcpyDeviceToHost);
dealError(err);
printFloatVec(h_dst, col, row);

native version

不利用share-memory，直接访global memory进行计算。这个版本比较好理解，不进行解释。

代码：

__global__ void transpose(float* src, float* dst, int row, int col){
int Idx = blockDim.x *blockIdx.x + threadIdx.x;
int Idy = blockDim.y *blockIdx.y + threadIdx.y;
//(Idx, Idy) =>(1D coordinate), col priority.
int indexIn = Idx + Idy * row;
int indexOut = Idy + Idx * col;
if(Idx < row && Idy < col){
dst[indexOut] = src[indexIn];
}
}

utilize share-memory version

思路：每个sm处理一块数据，每个sm内的线程先将数据从global memory拷贝进share memory，然后在share memory中完成转置，接着传回主存。这里，难点是地址的计算。

如何理解呢？我也没有想到好的解释。直接上代码吧。建议画图理解。

TODO

矩阵转置相关的论文（拓展）

分类：

方阵的转置，矩形矩阵的转置。

out-of-place， in-place转置。