在CUDA中具有共享MEM的非方面矩阵转置

Non-square matrix transpose with shared mem in CUDA

本文关键字:方面矩 转置 MEM 共享 CUDA      更新时间:2023-10-16

i试图为各种尺寸的cuda矩阵转置样品的变体变体。简而言之,我必须采用一个输入数组(double *a),然后将其写入一个较大矩阵(double *tab)的两个不同部分(您会注意到不同的偏移)。我正在以行 - 莫乔尔格式存储数据,因此我正在使用此宏来索引:

#define IDX2L(i,j,ld) (((i)*ld))+(j)) // 0 based index +row-major format

这是我使用的简单代码。

__global__ void cuda_a_Coalesced(double *tab, int tab_rows, int a_rows, double *a)
{
    __shared__  double tile[16*(16+1)]; 
    int col = threadIdx.x + blockIdx.x * blockDim.x;
    int row = threadIdx.y + blockIdx.y * blockDim.y;
    int col_2, row_2;
    int a_cols=tab_rows-a_rows; // tab_rows-a_rows is the number of columns of a
    int tab_cols=2*tab_rows+2;  // 2*tab_rows+2 is the number of columns of tab
    if( (col<a_cols) && (row<a_rows) ) 
    {
        // Load the data into shared mem
        tile[threadIdx.x+threadIdx.y*(16+1)]=a[IDX2L(row,col,a_cols)];
        // Normal copy (+ offsets)
        tab[IDX2L(row,col+tab_rows+a_rows,tab_cols)]= tile[threadIdx.x+threadIdx.y*(16+1)];
        // New idx
        col_2 = blockIdx.y * blockDim.y + threadIdx.x;
        row_2 = blockIdx.x * blockDim.x + threadIdx.y;
    }
    __syncthreads();
    if( (row_2<a_cols) && (col_2<a_rows) )
        // Transpose (+ other offsets)
        tab[IDX2L(row_2+a_rows,col_2+tab_rows,tab_cols)]= -tile[threadIdx.y+threadIdx.x*(16+1)];
}

启动参数是以下内容:

b1=(int)ceil((float)a_cols/16);
b2=(int)ceil((float)a_rows/16);
dim bck(b1,b2):dim th(16,16);
cuda_a_Coalesced<<<bck,th>>>(tab,tab_rows,a_rows,a);

无论大小如何,普通副本总是可以很好地执行。转置副本仅适用于块大小的整数倍数(如CUDA样本中)。当转置副本失败时,操作的某些部分是正确的,而其他部分则不正确,而我无法完全预测或跟踪。注意,因为想法是要在共享内存中更改索引,以便可以在输出矩阵中以合并形式写入旋转(由于行主要形式)。

有人可以告诉我代码仅适用于这种尺寸的原因?

有什么技巧可以解决这种情况?

问题是由于某些未定义的线程引起的,因为col_2row_2的值正在if()语句中分配,该语句中没有所有线程都在访问。

为了解决这种情况,当我们声明这些变量并删除已放置在上述if()中的同质计算时,我们可以给出col_2row_2的值

__shared__  double tile[16*(16+1)];
int col = threadIdx.x + blockIdx.x * blockDim.x;
int row = threadIdx.y + blockIdx.y * blockDim.y;
int col_2 = blockIdx.y * blockDim.y + threadIdx.x;
int row_2 = blockIdx.x * blockDim.x + threadIdx.y;
int a_cols=tab_rows-a_rows; 
int tab_cols=2*tab_rows+2;

因此,其余代码看起来像这样:

if( (col<a_cols) && (row<a_rows) ) 
{
    // Load the data into shared mem
    tile[threadIdx.x+threadIdx.y*(16+1)]=a[IDX2L(row,col,a_cols)];
    // Normal copy (+ offsets)
    tab[IDX2L(row,col+tab_rows+a_rows,tab_cols)]= tile[threadIdx.x+threadIdx.y*(16+1)];
}
__syncthreads();
if( (row_2<a_cols) && (col_2<a_rows) )
    // Transpose (+ other offsets)
    tab[IDX2L(row_2+a_rows,col_2+tab_rows,tab_cols)]= -tile[threadIdx.y+threadIdx.x*(16+1)];