备注

这个笔记是在自学过程的一些记录，中间包含了很多个人的理解以及问AI的回答，因为找不到很系统的教程，官方文档有时候不会讲的特别细，比方说某个指标的具体计算过程网上也没找到什么资料，所以在这里记录下，后续自己有了新的理解或者发现写的不对的地方都会在这不断更新，如有不对的地方，还请指正。

矩阵转置CPU朴素版本

这个没啥好说的，就是把A转置成B，A_ij=B_ji

void matrix_transponse_naive(const float* input, float* output, int rows, int cols) {
    for (int i = 0; i < rows * cols; i++) {
        output[i] = 0.0f;
    }

    for (int i = 0; i < rows; i++) {
        for (int j = 0; j < cols; j++) {
            output[j * rows + i] = input[i * cols + j];
        }
    }
}

矩阵转置CUDA朴素版本

CUDA朴素版本和CPU版本思路一致，只是分配了多个线程，每个线程负责一个元素的读取和写入，这里有个问题是我们在读取矩阵元素的时候是按行优先读取的，写入的时候是按列优先写入，导致读取是合并访存，但写入不是。

__global__ void matrix_transponse_naive_kernel(const float* input, float* output, int rows, int cols) {
    int tx = threadIdx.x;
    int ty = threadIdx.y;
    int bx = blockIdx.x;
    int by = blockIdx.y;

    int row_input = by * blockDim.y + ty;
    int col_input = bx * blockDim.x + tx;

    if (row_input < rows && col_input < cols) {
        int tid_input = row_input * cols + col_input;
        int tid_output = col_input * rows + row_input;
        output[tid_output] = input[tid_input];
    }

}

矩阵转置CUDA共享内存版本1

这里我们先把数据搬运到共享内存，但是非合并访存的问题依然存在

#define BLOCKDIM 16
__global__ void matrix_transponse_shared_kernel(const float* input, float* output, int rows, int cols) {
    __shared__ float smem[BLOCKDIM][BLOCKDIM];

    int tx = threadIdx.x;
    int ty = threadIdx.y;
    int bx = blockIdx.x;
    int by = blockIdx.y;

    int row = by * BLOCKDIM + ty;
    int col = bx * BLOCKDIM + tx;

    if (row < rows && col < cols) {
        smem[ty][tx] = input[row * cols + col];
    }
    __syncthreads();

    if (row < rows && col < cols) {
        output[col * rows + row] = smem[ty][tx];
    }

}

矩阵转置CUDA共享内存版本2

我们把上面的kernel做一个改进，在从全局内存读取数据写入到共享内存的时候，全局内存索引按行优先，然后从共享内存写入到全局内存的时候，全局内存索引也按行优先，相当于我把非合并访存的部分放在共享内存去做，由于共享内存属于片上内存，访问速度非常快，所以整体上耗时会进一步降低，
在全局内存中非合并：排队等候的是显存（DRAM/HBM），一次排队延迟是 200 ~ 400 个时钟周期。在共享内存中非合并（Bank Conflict）：排队发生在片上 SRAM，一次排队延迟通常只有几到十几个时钟周期。

__global__ void matrix_transponse_shared2_kernel(const float* input, float* output, int rows, int cols) {
    __shared__ float smem[BLOCKDIM][BLOCKDIM];

    int tx = threadIdx.x;
    int ty = threadIdx.y;
    int bx = blockIdx.x;
    int by = blockIdx.y;

    int row = by * BLOCKDIM + ty;
    int col = bx * BLOCKDIM + tx;

    if (row < rows && col < cols) {
        smem[ty][tx] = input[row * cols + col];
    }
    __syncthreads();

    int out_y = threadIdx.y + BLOCKDIM * blockIdx.x;
    int out_x = threadIdx.x + BLOCKDIM * blockIdx.y;
    int out_cols = rows;
    int out_rows = cols;
    if (out_x < out_cols && out_y < out_rows) {
        output[out_x + out_cols * out_y] =
            smem[threadIdx.x][threadIdx.y];
    }

}

但是这样带来了一个问题就是会产生Bank Conflicts，我们看下冲突是怎么产生的，以矩阵维度256256，blocksize=1616为例：

从上表我们可以看到一个warp中32个线程访问了4个bank，每个bank被8个线程同时访问，也就是8-way Bank Conflict，所以一个warp有8个wavefronts，理论上只需要一个。总共有256256/32=2048个warp，所以理论wavefronts=2048，实际wavefronts=20488=16384，Bank Conflicts=16384 - 2048 = 14336。

矩阵转置CUDA共享内存版本3

要消除Bank Conflicts通常可以使用Padding，也就是分配共享内存时额外+1

__global__ void matrix_transponse_shared3_kernel(const float* input, float* output, int rows, int cols) {
    __shared__ float smem[BLOCKDIM][BLOCKDIM + 1];

    int tx = threadIdx.x;
    int ty = threadIdx.y;
    int bx = blockIdx.x;
    int by = blockIdx.y;

    int row = by * BLOCKDIM + ty;
    int col = bx * BLOCKDIM + tx;

    if (row < rows && col < cols) {
        smem[ty][tx] = input[row * cols + col];
    }
    __syncthreads();

    int out_y = threadIdx.y + BLOCKDIM * blockIdx.x;
    int out_x = threadIdx.x + BLOCKDIM * blockIdx.y;
    int out_cols = rows;
    int out_rows = cols;
    if (out_x < out_cols && out_y < out_rows) {
        output[out_x + out_cols * out_y] =
            smem[threadIdx.x][threadIdx.y];
    }

}

我们分析下这段代码的bank，从下表可以看到共享内存空间变为1617后，一个warp中32个线程访问了31个线程，只有bank0有两个线程访问，所以一个warp也就只有一个Bank Conflict。

总共有256256/32=2048个warp，所以理论wavefronts=2048，实际wavefronts=2048*2=4096，Bank Conflicts=4096 - 2048 = 2048。

不同kernel耗时

在矩阵维度为2560025600个Float、blocksize=1616下对比不同kernel耗时结果如下