我们如何访问 CUDA 中的 3D 数组的列？答案

【问题标题】：How do we access the column of a 3D array in CUDA?我们如何访问 CUDA 中的 3D 数组的列？
【发布时间】：2020-05-12 00:53:25
【问题描述】：

    mod=SourceModule("""

   __global__ void mat_ops(float *A,float *B)
  {   /*formula to get unique thread index*/
      int thrd= blockIdx.x*blockDim.x*blockDim.y+threadIdx.y*blockDim.x+threadIdx.x;
      B[]=A[];
   }    
   """)
        func = mod.get_function("mat_ops")
        func(A_k, B_k, grid=(3,1,1),block=(4,4,1))

我有两个 3D 数组 float *A 和 float *B，在这个 PyCUDA 内核中每个大小为 4 X 4 X 3。我在这里要做的是，逐列而不是逐行遍历 3D 数组。我正在使用 2D 块的 1D 网格。我该怎么做呢？

【问题讨论】：

在不知道数组顺序的情况下无法回答这个问题，你没有描述过
数组按行主要顺序排列。您也可以在此处使用单词矩阵作为数组。换句话说，我所说的数组，我的意思是这里的矩阵
在数学上，将 3D 数组称为矩阵是错误的，因为它不是一个。所以你想在这种情况下遍历numpy数组的中轴，还是第一个？
我想访问 CUDA 内核中 3D 数组的中轴。另外，感谢您在术语上纠正我。

标签： python numpy cuda pycuda

【解决方案1】：

为此，您需要向 CUDA 内核描述内存中数组的布局，并且您需要使用主机端提供的步幅在内核中进行正确的索引计算。一个简单的方法是在 CUDA 中定义一个小的帮助类，它隐藏了大部分索引并提供了一个简单的索引语法。例如：

from pycuda import driver, gpuarray
from pycuda.compiler import SourceModule
import pycuda.autoinit
import numpy as np

mod=SourceModule("""

   struct stride3D
   {
       float* p;
       int s0, s1;

       __device__
       stride3D(float* _p, int _s0, int _s1) : p(_p), s0(_s0), s1(_s1) {};

       __device__
       float operator  () (int x, int y, int z) const { return p[x*s0 + y*s1 + z]; };

       __device__
       float& operator () (int x, int y, int z) { return p[x*s0 + y*s1 + z]; };
   };

   __global__ void mat_ops(float *A, int sA0, int sA1, float *B, int sB0, int sB1)
   {
       stride3D A3D(A, sA0, sA1);
       stride3D B3D(B, sB0, sB1);

       int xidx = blockIdx.x;
       int yidx = threadIdx.x;
       int zidx = threadIdx.y;

       B3D(xidx, yidx, zidx) = A3D(xidx, yidx, zidx);
   }    
   """)

A = 1 + np.arange(0, 4*4*3, dtype=np.float32).reshape(4,4,3)
B = np.zeros((5,5,5), dtype=np.float32)
A_k = gpuarray.to_gpu(A)
B_k = gpuarray.to_gpu(B)

astrides = np.array(A.strides, dtype=np.int32) // A.itemsize
bstrides = np.array(B.strides, dtype=np.int32) // B.itemsize

func = mod.get_function("mat_ops")
func(A_k, astrides[0], astrides[1], B_k, bstrides[0], bstrides[1], grid=(4,1,1),block=(4,3,1))
print(B_k[:4,:4,:3])

这里我选择使源数组和目标数组大小不同，只是为了表明代码是通用的，只要块大小足够，它就可以适用于任何大小的数组。请注意，这里没有在设备代码方面检查数组边界，您需要为非平凡的示例添加它。

~~还请注意，这对于 fortran 和 C 有序 numpy 数组都应该正确工作，因为它直接使用 numpy 步幅值。但是，由于内存合并问题，CUDA 方面的性能会受到影响。~~

注意：如果不扩展帮助程序类以在所有维度上跨步并更改内核以接受输入和输出数组的所有维度的跨步，这将不适用于 fortran 和 C 排序。从性能的角度来看，最好为 fortran 和 C 有序数组编写单独的辅助类。

【讨论】：