python - PyCUDA 核函数

Question

我是 PyCUDA 的新手，正在浏览 PyCUDA 网站上的一些示例。我正在尝试弄清楚某些代码行背后的逻辑，如果有人解释了它背后的想法，我将不胜感激。

以下代码片段来自 PyCUDA 网站。函数定义里面，没看懂

int idx = threadIdx.x + threadIdx.y*4;

如何使用上面的行来计算数组的索引。为什么threadIdx.x和threadIdx.y加在一起，为什么threadIdx.y乘以4。

对于 GPU 的函数调用，为什么 block 定义为 5,5,1。因为它是一个 5x5 元素的数组，所以根据我的理解， block 大小应该是 5,5 而不是 5,5,1。

import pycuda.driver as cuda
import pycuda.autoinit
from pycuda.compiler import SourceModule
import numpy
a = numpy.random.randn(5,5)
a = a.astype(numpy.float32)
a_gpu = cuda.mem_alloc(a.nbytes)
cuda.memcpy_htod(a_gpu, a)
mod = SourceModule("""
__global__ void doubleMatrix(float *a)
{
int idx = threadIdx.x + threadIdx.y*4;
a[idx] *= 2;
}
""")
func = mod.get_function("doubleMatrix")
func(a_gpu, block=(5,5,1))
a_doubled = numpy.empty_like(a)
cuda.memcpy_dtoh(a_doubled, a_gpu)
print ("ORIGINAL MATRIX")
print a
print ("DOUBLED MATRIX AFTER PyCUDA EXECUTION")
print a_doubled

Answer 1

您发布的示例似乎来自(或抄袭)一本名为“Python Parallel Programming Cookbook”的书，直到五分钟前我才听说过这本书。老实说，如果我是那本书的作者，我会为包含这样一个骇人听闻的错误示例而感到羞愧。

这是对您发布的内容及其输出的一个小修改:

import pycuda.driver as cuda
import pycuda.autoinit
from pycuda.compiler import SourceModule
import numpy
a = numpy.random.randn(5,5)
a = a.astype(numpy.float32)
a_gpu = cuda.mem_alloc(a.nbytes)
cuda.memcpy_htod(a_gpu, a)
mod = SourceModule("""
__global__ void doubleMatrix(float *a)
{
     int idx = threadIdx.x + threadIdx.y*4;
     a[idx] *= 2.f;
}
""")
func = mod.get_function("doubleMatrix")
func(a_gpu, block=(5,5,1))
a_doubled = numpy.empty_like(a)
cuda.memcpy_dtoh(a_doubled, a_gpu)
print a_doubled - 2.0*a

[警告:Python 2 语法]

In [2]: %run matdouble.py
[[ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.         -0.62060976  0.49836278 -1.60820103  1.71903515]]

即代码没有按预期工作，这可能是您感到困惑的根源。

在 this very recent answer 中描述了寻址存储在线性内存中的多维数组(如 numpy 数组)的正确方法。 .任何明智的程序员都会像这样在您的示例中编写内核:

__global__ void doubleMatrix(float *a, int lda)
{
     int idx = threadIdx.x + threadIdx.y * lda;
     a[idx] *= 2.f;
}

以便将数组的前导维度作为参数传递给内核(在本例中应为 5，而不是 4)。这样做会得到以下结果:

import pycuda.driver as cuda
import pycuda.autoinit
from pycuda.compiler import SourceModule
import numpy
a = numpy.random.randn(5,5)
a = a.astype(numpy.float32)
a_gpu = cuda.mem_alloc(a.nbytes)
cuda.memcpy_htod(a_gpu, a)
mod = SourceModule("""
__global__ void doubleMatrix(float *a, int lda)
{
     int idx = threadIdx.x + threadIdx.y * lda;
     a[idx] *= 2.f;
}
""")
func = mod.get_function("doubleMatrix")
lda = numpy.int32(a.shape[-1])
func(a_gpu, lda, block=(5,5,1))
a_doubled = numpy.empty_like(a)
cuda.memcpy_dtoh(a_doubled, a_gpu)
print a_doubled - 2.0*a

产生预期结果:

In [3]: %run matdouble.py
[[ 0.  0.  0.  0.  0.]
 [ 0.  0.  0.  0.  0.]
 [ 0.  0.  0.  0.  0.]
 [ 0.  0.  0.  0.  0.]
 [ 0.  0.  0.  0.  0.]]

For the function call to the GPU why is the block defined as 5,5,1. Since it is an array of 5x5 elements so in my understanding the block size should be 5,5 instead of 5,5,1.

在 CUDA 中，所有 block 都隐式具有三个维度。 (5,5) 的 block 大小与 (5,5,1) 的 block 大小相同。最后一个维度可以忽略，因为它是一个维度(即 block 中的所有线程都有 threadIdx.z = 1)。您不应该陷入的陷阱是将 CUDA block 或网格的维度与输入数组的维度混为一谈。有时让它们相同很方便，但同样没有必要甚至不建议这样做。对于此示例(假设行主要存储顺序)，以 BLAS 样式正确编写的内核可能如下所示:

__global__ void doubleMatrix(float *a, int m, int n, int lda)
{
     int col = threadIdx.x + blockIdx.x * blockDim.x;
     int row = threadIdx.y + blockDim.y * blockDim.y;

     for(; row < m; row += blockDim.y * gridDim.y) {
         for(; col < n; col += blockDim.x * gridDim.x) {
             int idx = col + row * lda;
             a[idx] *= 2.f;
         }
    }
}

[注:在浏览器中编写，未经编译或测试]

这里 任何 合法的 block 和网格维度将正确处理任何大小的输入数组元素总数将适合带符号的 32 位整数。如果你运行太多线程，有些线程什么都不做。如果运行的线程太少，一些线程将处理多个数组元素。如果您运行一个与输入数组具有相同维度的网格，每个线程将只处理一个输入，正如您正在研究的示例中的意图。如果您想阅读有关如何选择最合适的 block 和网格大小的信息，我建议您开始 here .