0

我刚刚开始 CUDA 编程，并且正在学习有关矩阵乘法的内核设计的知识。我复制了网上找到的主要代码，并添加了我的部分来实现 A(MxM) 和 B(MxN) 的矩阵乘法

#include <iostream>
#include <iomanip>
#include <fstream>
#include <vector>
#include <cuda_runtime.h>
#include <cuComplex.h>
#include <cusp/complex.h>

using namespace std;

const int M=55, N=73;

typedef cusp::complex<double> Complex;

__global__ void kernelFunc(Complex* ad, Complex* bd, Complex* cd, int n) 
{
  int x = (blockIdx.x * blockDim.x) + threadIdx.x;
  int y = (blockIdx.y * blockDim.y) + threadIdx.y;

  if(x < n && y < n) 
  {   
    Complex v = Complex(0.0, 0.0);
    for(int i=0; i<n; i++) v += ad[y * n + i] * bd[i * n + x];
    cd[y * n + x] = v;
  }
}

int main(int argc, char *argv[])
{
  std::vector< Complex > _A(M*M);
  std::vector< Complex > _B(M*N);
  Complex *A, *B, *C;

  cudaMalloc((void**)&A, M*M*sizeof(Complex));
  cudaMalloc((void**)&B, M*N*sizeof(Complex));
  cudaMalloc((void**)&C, M*N*sizeof(Complex));

  for (int i=0; i<M*M; i++) _A[i] = Complex((double)i, (double)i);
  for (int i=0; i<M*N; i++) _B[i] = Complex(1.0, 0.0);

  cudaMemcpy( A, &_A[0], (M*M)*sizeof(Complex), cudaMemcpyHostToDevice );
  cudaMemcpy( B, &_B[0], (M*N)*sizeof(Complex), cudaMemcpyHostToDevice );

  dim3 block(32, 32);           
  dim3 grid((M+31)/32, (M+31)/32);

  kernelFunc<<<grid, block>>>(A, B, C, M);
  cudaMemcpy(&_B[0], &C[0], (M*N)*sizeof(Complex), cudaMemcpyDeviceToHost);

  cudaFree(A);
  cudaFree(B);
  cudaFree(C);

  return 0;
}

但是网上说涉及的矩阵一定是方阵，是不是说这段代码不能用在任意维数的矩阵上呢？我不明白如何定义块数和网格数以适应我的问题。我的问题中的矩阵的维度为 MxM，其中 M 为奇数。我尝试上面的小矩阵代码，它似乎有效，但我必须将它应用于相当大的矩阵。我不知道这是否也适用于大矩阵。

4

1 回答 1

2

给定的代码仅适用于方阵。需要稍作修改才能使其通用化。内核可以修改如下：

__global__ void kernelFunc(Complex* ad, Complex* bd, Complex* cd, int m1, int n1, int n2) 
{
  int x = (blockIdx.x * blockDim.x) + threadIdx.x;
  int y = (blockIdx.y * blockDim.y) + threadIdx.y;

  if(x < n2 && y < m1) 
  {   
    Complex v = Complex(0.0, 0.0);
    for(int i=0; i<n1; i++) v += ad[y * n1 + i] * bd[i * n2 + x];
    cd[y * n2 + x] = v;
  }
}

其中m1是第一个矩阵的行，是第一个矩阵n1的列，n2是第二个矩阵的列。网格大小将相应修改为：

dim3 grid((N+31)/32, (M+31)/32);

最后，内核调用：

kernelFunc<<<grid, block>>>(A, B, C, M, M, N);
于 2013-09-11T18:49:05.680 回答