如何推广方阵乘法以处理任意维度
how to generalize square matrix multiplication to handle arbitrary dimensions
我已经编写了这个程序,但我在理解如何通过在内核调用行中使用 dim3 变量来使用多个块时遇到了一些麻烦。当我进行 1000*1000 矩阵乘法时,这段代码工作正常,但没有得到较低维度的正确答案,如 100*100 、 200*200.
#include <stdio.h>
#include <cuda.h>
#define width 1000
__global__ void kernel(int *a,int *b,int *c)
{
int tx = threadIdx.x + blockIdx.x*blockDim.x;
int ty = threadIdx.y + blockIdx.y*blockDim.y;
int sum=0,k;
for(k=0;k<(width);++k)
{
sum += a[ty*width +k]*b[k*width + tx];
}
c[ty*width + tx] = sum;
}
int main()
{
int a[width*width],c[width*width],b[width*width];
int *dev_a,*dev_b,*dev_c;
int i,count=0;
int size = (width*width)*sizeof(int);
for(i=0;i<(width*width);i++)
{
a[i] = 1;
b[i] = 1;
}
cudaMalloc((void **)&dev_a,size);
cudaMalloc((void **)&dev_b,size);
cudaMalloc((void **)&dev_c,size);
cudaMemcpy(dev_a,&a,size,cudaMemcpyHostToDevice);
cudaMemcpy(dev_b,&b,size,cudaMemcpyHostToDevice);
dim3 dimBlock(20,20);
dim3 blockID(50,50);
kernel<<<blockID,dimBlock>>>(dev_a,dev_b,dev_c);
cudaMemcpy(&c,dev_c,size,cudaMemcpyDeviceToHost);
for(i=0;i<(width*width);i++)
{
count++;
if(count == (width+1))
{
count = 1;
printf("\n");
}
printf("%d ",c[i]);
}
printf("\n");
return 0;
}
此代码适用于非常具体的维度,但不适用于其他维度。
当 width
正好等于您的块维度(线程数 - 您显示的代码中的 20)和您的网格维度(块数 -您显示的代码中有 50 个)。
因此,当 width
为 20*50 (1000) 时,它将如图所示工作。但是,如果我将 width
更改为其他值(例如 800)并且不进行任何其他更改,您的代码将无法运行。但是,在 800 的情况下,我可以通过将网格维度从 50 更改为 40 来让您的代码正常工作,然后 width
= 800 = 20 *40.
但是如果我需要将两个 width
799 的矩阵相乘怎么办?我无法想出与 width
完全匹配的网格和块尺寸的乘积。
这是 CUDA 编程中的一个相当标准的问题 - 我无法想出方便的块和网格尺寸来完全匹配我的工作(即数据)大小,如果我启动太多(threads/blocks)东西好像不行。
要解决这个问题,我们必须做两件事:
- 确保至少启动足够多的线程(线程块)来覆盖整个数据集
- 在内核中添加条件代码,以便只有对应于有效数据的线程才做真正的工作。
为了解决上面的第 1 项,我们将网格维度计算修改为如下内容:
dim3 dimBlock(16,16);
dim3 blockID((width+dimBlock.x-1)/dimBlock.x,(width+dimBlock.y-1)/dimBlock.y);
为了解决上面的第 2 项,我们修改了内核代码以根据线程是否对应于有效数据来调节线程行为:
__global__ void kernel(int *a,int *b,int *c, int mwidth)
{
int tx = threadIdx.x + blockIdx.x*blockDim.x;
int ty = threadIdx.y + blockIdx.y*blockDim.y;
if ((tx<mwidth)&&(ty<mwidth)){
int sum=0,k;
for(k=0;k<(mwidth);++k)
{
sum += a[ty*mwidth +k]*b[k*mwidth + tx];
}
c[ty*mwidth + tx] = sum;}
}
并且因为我们用新参数修改了内核,所以我们必须在调用时传递该参数:
kernel<<<blockID,dimBlock>>>(dev_a,dev_b,dev_c, width);
这应该是逻辑扩展您显示的处理 "arbitrary" 维度的代码所需要的。我还建议在您遇到 CUDA 代码问题时添加 proper cuda error checking。
我已经编写了这个程序,但我在理解如何通过在内核调用行中使用 dim3 变量来使用多个块时遇到了一些麻烦。当我进行 1000*1000 矩阵乘法时,这段代码工作正常,但没有得到较低维度的正确答案,如 100*100 、 200*200.
#include <stdio.h>
#include <cuda.h>
#define width 1000
__global__ void kernel(int *a,int *b,int *c)
{
int tx = threadIdx.x + blockIdx.x*blockDim.x;
int ty = threadIdx.y + blockIdx.y*blockDim.y;
int sum=0,k;
for(k=0;k<(width);++k)
{
sum += a[ty*width +k]*b[k*width + tx];
}
c[ty*width + tx] = sum;
}
int main()
{
int a[width*width],c[width*width],b[width*width];
int *dev_a,*dev_b,*dev_c;
int i,count=0;
int size = (width*width)*sizeof(int);
for(i=0;i<(width*width);i++)
{
a[i] = 1;
b[i] = 1;
}
cudaMalloc((void **)&dev_a,size);
cudaMalloc((void **)&dev_b,size);
cudaMalloc((void **)&dev_c,size);
cudaMemcpy(dev_a,&a,size,cudaMemcpyHostToDevice);
cudaMemcpy(dev_b,&b,size,cudaMemcpyHostToDevice);
dim3 dimBlock(20,20);
dim3 blockID(50,50);
kernel<<<blockID,dimBlock>>>(dev_a,dev_b,dev_c);
cudaMemcpy(&c,dev_c,size,cudaMemcpyDeviceToHost);
for(i=0;i<(width*width);i++)
{
count++;
if(count == (width+1))
{
count = 1;
printf("\n");
}
printf("%d ",c[i]);
}
printf("\n");
return 0;
}
此代码适用于非常具体的维度,但不适用于其他维度。
当 width
正好等于您的块维度(线程数 - 您显示的代码中的 20)和您的网格维度(块数 -您显示的代码中有 50 个)。
因此,当 width
为 20*50 (1000) 时,它将如图所示工作。但是,如果我将 width
更改为其他值(例如 800)并且不进行任何其他更改,您的代码将无法运行。但是,在 800 的情况下,我可以通过将网格维度从 50 更改为 40 来让您的代码正常工作,然后 width
= 800 = 20 *40.
但是如果我需要将两个 width
799 的矩阵相乘怎么办?我无法想出与 width
完全匹配的网格和块尺寸的乘积。
这是 CUDA 编程中的一个相当标准的问题 - 我无法想出方便的块和网格尺寸来完全匹配我的工作(即数据)大小,如果我启动太多(threads/blocks)东西好像不行。
要解决这个问题,我们必须做两件事:
- 确保至少启动足够多的线程(线程块)来覆盖整个数据集
- 在内核中添加条件代码,以便只有对应于有效数据的线程才做真正的工作。
为了解决上面的第 1 项,我们将网格维度计算修改为如下内容:
dim3 dimBlock(16,16);
dim3 blockID((width+dimBlock.x-1)/dimBlock.x,(width+dimBlock.y-1)/dimBlock.y);
为了解决上面的第 2 项,我们修改了内核代码以根据线程是否对应于有效数据来调节线程行为:
__global__ void kernel(int *a,int *b,int *c, int mwidth)
{
int tx = threadIdx.x + blockIdx.x*blockDim.x;
int ty = threadIdx.y + blockIdx.y*blockDim.y;
if ((tx<mwidth)&&(ty<mwidth)){
int sum=0,k;
for(k=0;k<(mwidth);++k)
{
sum += a[ty*mwidth +k]*b[k*mwidth + tx];
}
c[ty*mwidth + tx] = sum;}
}
并且因为我们用新参数修改了内核,所以我们必须在调用时传递该参数:
kernel<<<blockID,dimBlock>>>(dev_a,dev_b,dev_c, width);
这应该是逻辑扩展您显示的处理 "arbitrary" 维度的代码所需要的。我还建议在您遇到 CUDA 代码问题时添加 proper cuda error checking。