如何为这个函数编写cuda内核？

Question

我有一个 C 内核，我想将其转换为 cuda 内核。然而，我在为这个函数编写内核时遇到了麻烦，因为它由多个嵌套循环组成。我的方法是将这个函数分成三个不同的内核，每个操作一个。 这是原来的功能：

void  kernel(int N, double B[N][N]) {

  for (int i = 0; i < N; i++) {
    for (int j = 0; j < i; j++) {
      for (int k = 0; k < j; k++) {
          B[i][j] = B[i][j] - (B[i][k] * B[k][j]);
      }
       B[i][j] = B[i][j] / B[j][j];
    }
    for (int j = i; j < N; j++) {
      for (int k = 0; k < i; k++) {
         B[i][j] = B[i][j] - (B[i][k] * B[k][j]);
      }
    }
  }
}

这是我尝试将其转换为 3 个不同的内核：

__global__ void kernal_1(int N, double* B){

 int i = blockIdx.x * blockDim.x + threadIdx.x;
 int j = blockIdx.y * blockDim.y + threadIdx.y;
 int k;

    if (i < N && j<i) {
        for (k = 0; k < j; k++) {
             B[i][j] = B[i][j] - (B[i][k] * B[k][j]);
            } 
    }

}

__global__ void  kernal_2(int N, double* B){

 int i = blockIdx.x * blockDim.x + threadIdx.x;
 int j = blockIdx.y * blockDim.y + threadIdx.y;

    if (i < N && j<i) {
      B[i][j] = B[i][j] / B[j][j];

    }

}


__global__ void kernal_3(int N, double* B){

int i = blockIdx.x * blockDim.x + threadIdx.x;
int j = blockIdx.y * blockDim.y + threadIdx.y;
int k;

    if (i < N && j<=i) {
    for (k = 0; k < i; k++) {
         B[i][j] = B[i][j] - (B[i][k] * B[k][j]);
        }    
       
}}

我知道我的尝试是错误的，因为我没有得到正确的结果。有什么意见或建议吗？顺便说一句，这是用于 LU 分解的。我见过人们采用不同的方式，其中外部 i 循环没有在内核中实现，而是在调用内核的主函数中使用，内核在每种情况下都有 2 个循环

编辑：包含内存访问模式

Answer 1

第一个答案是，你可以尝试CUDA中的动态并行，你可以在google上的CUDA编程指南中搜索相关内容。 （但对于 CUDA 初学者来说可能太难了）

第二个想法是，不需要将它们分成三个内核，只需使用一个内核，因为不同的内核将数据从全局传输到寄存器，再传输到全局再到寄存器，这绝对是非常非常慢的。使用CUDA就像一层for，所以你只需要考虑如何将任务分配给所有线程即可：

    for (int j = 0; j < i; j++) {
      for (int k = 0; k < j; k++) {
         // operation 1 
      }
      //operation 2
    }
    for (int j = i; j < N; j++) {
      for (int k = 0; k < i; k++) {
        // operation 3 
      }
    }

如果第一个迭代器，“i = 0; i < N; i++" is large enough to fill the GPU, just use for loop(for j and k) within the kernel, it is OK. What you should do is fully use GPU's capability. How to know whether N is large enough? So you firstly write a kernel, and use this 占用计算器（使用 nsight 计算查找共享内存和寄存器值）找出扭曲占用率，并获取一个 SM 中可以活动的块数，并乘以你拥有的SM总数。然后你就得到了。