2016-04-13 149 views
2

我已经写了下面的代码来将两个4x4矩阵求和。Cuda矩阵加法

#include<stdio.h> 
#include<stdlib.h> 
#include<math.h> 

__global__ void Matrix_add(double* a, double* b, double* c,int n) 
{ 
    int row = blockIdx.x * blockDim.x + threadIdx.x; 
    int col = blockIdx.y * blockDim.y + threadIdx.y; 
    int index = row * n + col; 
    if(col<n && row <n) 
     c[index] = a[index] + b[index]; 
} 
int main() 
{ 

int n=4; 

double **h_a; 
double **h_b; 
double **h_c; 
double *d_a, *d_b, *d_c; 

int size = n*n*sizeof(double); 

h_a = (double **) malloc(n*sizeof(double*)); 
h_b = (double **) malloc(n*sizeof(double*)); 
h_c = (double **) malloc(n*sizeof(double*)); 

cudaMalloc((void**)&d_a,size); 
cudaMalloc((void**)&d_b,size); 
cudaMalloc((void**)&d_c,size); 

int t=0; 
for (t=0;t<n;t++) 
{ 
     h_a[t]= (double *)malloc(n*sizeof(double)); 
     h_b[t]= (double *)malloc(n*sizeof(double)); 
     h_c[t]= (double *)malloc(n*sizeof(double)); 
} 

int i=0,j=0; 

for(i=0;i<n;i++) 
{ 
for(j=0;j<n;j++) 
    { 
     h_a[i][j]=sin(i)*sin(i); 
     h_b[i][j]=cos(i)*cos(i); 
    } 
} 

cudaMemcpy(d_a,h_a+n,size,cudaMemcpyHostToDevice); 
cudaMemcpy(d_b,h_b+n,size,cudaMemcpyHostToDevice); 

dim3 dimBlock(4,4); 
dim3 dimGrid(1,1); 
Matrix_add<<<dimGrid, dimBlock>>>(d_a,d_b,d_c,n); 
cudaMemcpy(h_c+n,d_c,size,cudaMemcpyDeviceToHost); 

for(i=0;i<n;i++) 
{ 
    for(j=0;j<n;j++) 
     { 
      printf("%f",h_c[i][j]); 
      printf("\t"); 
     } 
    printf("\n"); 
} 

for(i=0;i<n;i++) 
    { 
     free(h_a[i]); 
     free(h_b[i]); 
     free(h_c[i]); 
} 
free(h_a); 
free(h_b); 
free(h_c); 
cudaFree(d_a); cudaFree(d_b); cudaFree(d_c); 
return 0; 

} 

结果该相加应该是一个2×2全1矩阵,但在结果矩阵的所有元素都为0。此外,我得到这个消息得到结果后:

分段故障(核心倾销)

任何人都可以请帮我找出问题。

谢谢

回答

5

您的主机阵列(h_a,数组h_b​​,h_c)不在内存中连续的,所以你的初始cudaMemcpy()调用将读取垃圾变成GPU内存(在你的情况显然是零)。

原因是你的主机数组实际上并不是平坦的,而是表示为指针数组。我猜在C中伪造二维数组?无论如何,您需要更加小心cudaMemcpy(),并逐行复制主机数组,或在主机上使用平面表示。

+0

很好的答案。正确的解决方案是使用平面表示。该单独分配的行数组很乱,容易出错并且效率低下。你会习惯'h [i * row_stride + j]'(或者将它隐藏在宏的后面),并且你的簿记代码将大大缩小。 – Peter