CUDA并行化依赖2D阵列

我有以下形式的示例循环。请注意，我的psi[i][j]依赖于psi[i+1][j], psi[i-1][j], psi[i][j+1] and psi[i][j-1]，我只能计算内部矩阵的psi。现在，我尝试在CUDA中编写此代码，但结果与顺序不同。CUDA并行化依赖2D阵列

for(i=1;i<=leni-2;i++) 
for(j=1;j<=lenj-2;j++){ 
    psi[i][j]=(omega[i][j]*(dx*dx)*(dy*dy)+(psi[i+1][j]+psi[i-1][j])*(dy*dy)+(psi[i][j+1]+psi[i][j-1])*(dx*dx))/(2.0*(dx*dx)+2.0*(dy*dy)); 
}

这是我的CUDA格式。

//KERNEL 
__global__ void ComputePsi(double *psi, double *omega, int imax, int jmax) 
{ 
int x = blockIdx.x; 
int y = blockIdx.y; 
int i = (jmax*x) + y; 
double beta = 1; 
double dx=(double)30/(imax-1); 
double dy=(double)1/(jmax-1); 

if((i)%jmax!=0 && (i+1)%jmax!=0 && i>=jmax && i<imax*jmax-jmax){ 
    psi[i]=(omega[i]*(dx*dx)*(dy*dy)+(psi[i+jmax]+psi[i-jmax])*(dy*dy)+(psi[i+1]+psi[i-1])*(dx*dx))/(2.0*(dx*dx)+2.0*(dy*dy)); 
} 
} 


//Code 
cudaMalloc((void **) &dev_psi, leni*lenj*sizeof(double)); 
cudaMalloc((void **) &dev_omega, leni*lenj*sizeof(double)); 
cudaMemcpy(dev_psi, psi, leni*lenj*sizeof(double),cudaMemcpyHostToDevice); 
cudaMemcpy(dev_omega, omega, leni*lenj*sizeof(double),cudaMemcpyHostToDevice); 
dim3 grids(leni,lenj); 
for(iterpsi=0;iterpsi<30;iterpsi++)   
    ComputePsi<<<grids,1>>>(dev_psi, dev_omega, leni, lenj);

其中psi[leni][lenj] and omega[leni][lenj]和双数组。

问题在于顺序问题，CUDA代码给出了不同的结果。代码中是否需要修改？

来源

2014-04-07 unkn0wn

无关提示：我总是将输入缓冲区声明为const指针，所以我不会搞砸。 – texasbruce

相关提示：始终将输入和输出缓冲区分开，并且不要写入全局输入缓冲区 – texasbruce

您正在全局内存中工作，并且您正在更改psi条目，而其他线程可能需要旧值。只需将新迭代的值存储在单独的变量中。但请记住，每次迭代后必须交换变量！更复杂的方法是使用共享内存和空间域分配到独立线程的解决方案。只需谷歌的CUDA教程解决热/扩散方程，你会明白。

来源

2014-04-07 21:17:33 OnWhenReady

感谢您的回复，最初我在每次迭代后交换变量。但是由于我的迭代次数非常多，有人建议在cuda上交换需要很多时间，所以我正在寻找更快的替代方法。你的第二个解决方案看起来很有趣，我会看看它。谢谢。 – unkn0wn

你不必真的交换。只需在你的内核中添加一个bool，并在每次迭代之后交换你的主机上的bool。那么你的内核必须实现如 'if（swap）{temp_psi = ... psi ...;} else {psi = ... temp_psi ...}'。 – OnWhenReady

Helllo在那里，我试图寻找你第一个建议的交换方法，但可以找到任何解决方案值得我。您能否提出任何适合的解决方案或链接详细说明使用交换变量来解决gauss-seidel问题。谢谢。 – unkn0wn

for(i=1;i<=leni-2;i++) 
    for(j=1;j<=lenj-2;j++){ 
    psi[i][j]= (omega[i][j]*(dx*dx)*(dy*dy) + 
       (psi[i+1][j]+psi[i-1][j]) * (dy*dy) + 
       (psi[i][j+1]+psi[i][j-1]) * (dx*dx) 
       )/(2.0*(dx*dx)+2.0*(dy*dy)); 
}

我认为，这个内核是不正确的顺序之一：的psi[i][j]的值取决于操作的顺序放在这里 - 所以你将使用不更新psi[i+1][j]和psi[i][j+1]，但psi[i-1][j]和psi[i][j-1]已更新在这个扫描。

请确保使用CUDA的结果将会不同，其中操作顺序不同。

为了执行这样的排序，如果可能的话，您需要插入如此多的同步，这可能对CUDA来说不值得。它真的是你需要做的吗？

来源

2014-04-08 03:17:09 Sigismondo

使用这个内核的同步根本无济于事。正如我所建议的OP应该使用交换变量。这种核函数（思想）是许多人（包括我）使用CUDA解决的微分方程的典型例子。并行解决PDE/ODE是CUDA确实存在的原因;-) – OnWhenReady

是的，我认为这是错误的顺序 - 但是如果它真的是他打算做的......通常你需要编写后向有限差分PDE内核为'psi_new [i] [j] = f（psi）'，然后用'psi_new'交换'psi'。并行和顺序。让我们看看OP如何看待它。 – Sigismondo

这正是我说的话:-) – OnWhenReady

CUDA并行化依赖2D阵列

回答

相关问题