内核并行

我目前正在尝试使用GPU并行计算处理雷达信号，需要在200ms以内同时完成6批数据的一系列处理。每批的数据量是4096×512，我使用的显卡是Tesla P40.
我的gridDim是（32，16，16），blockDim是（8，8，4）
目前我遇到的问题是：有一个计算环节，我将数据排列成二维的矩阵，需要对每个数据单元，取其上下左右各16个数据，共64个数据求和取平均，计算量很大。为了处理边界问题，我还使用了纹理数组，但即使如此，计算时间还是超出预算。

该部分的代码如下：
__global__ void CFAR0(float *cfar, float *real)
{
   int bid = blockIdx.x + blockIdx.y*gridDim.x + blockIdx.z*gridDim.x*gridDim.y;
   int tid = bid * blockDim.x * blockDim.y*blockDim.z + threadIdx.x + threadIdx.y * blockDim.x + threadIdx.z*blockDim.x*blockDim.y;
   float x = tid % DistanceDim;
   float y = tid / DistanceDim;
   float average1 = (tex2D(texRef01, x, y - 9) + tex2D(texRef01, x, y - 10)
      + tex2D(texRef01, x, y - 11) + tex2D(texRef01, x, y - 12)
      + tex2D(texRef01, x, y - 13) + tex2D(texRef01, x, y - 14)
      + tex2D(texRef01, x, y - 15) + tex2D(texRef01, x, y - 16)
      + tex2D(texRef01, x, y - 17) + tex2D(texRef01, x, y - 18)
      + tex2D(texRef01, x, y - 19) + tex2D(texRef01, x, y - 20)
      + tex2D(texRef01, x, y - 21) + tex2D(texRef01, x, y - 22)
      + tex2D(texRef01, x, y - 23) + tex2D(texRef01, x, y - 24)) / 16;

   cfar[tid] = (real[tid] >= average1 * threshold) ? real[tid] : 0;
   if (cfar[tid] != 0){
      float average2 = (tex2D(texRef01, x, y + 9) + tex2D(texRef01, x, y + 10)
         + tex2D(texRef01, x, y + 11) + tex2D(texRef01, x, y + 12)
         + tex2D(texRef01, x, y + 13) + tex2D(texRef01, x, y + 14)
         + tex2D(texRef01, x, y + 15) + tex2D(texRef01, x, y + 16)
         + tex2D(texRef01, x, y + 17) + tex2D(texRef01, x, y + 18)
         + tex2D(texRef01, x, y + 19) + tex2D(texRef01, x, y + 20)
         + tex2D(texRef01, x, y + 21) + tex2D(texRef01, x, y + 22)
         + tex2D(texRef01, x, y + 23) + tex2D(texRef01, x, y + 24)) / 16;
      cfar[tid] = (cfar[tid] >= average2 * threshold) ? cfar[tid] : 0;
   }
   if (cfar[tid] != 0){
      float average3 = (tex2D(texRef01, x - 9, y) + tex2D(texRef01, x - 10, y)
         + tex2D(texRef01, x - 11, y) + tex2D(texRef01, x - 12, y)
         + tex2D(texRef01, x - 13, y) + tex2D(texRef01, x - 14, y)
         + tex2D(texRef01, x - 15, y) + tex2D(texRef01, x - 16, y)
         + tex2D(texRef01, x - 17, y) + tex2D(texRef01, x - 18, y)
         + tex2D(texRef01, x - 19, y) + tex2D(texRef01, x - 20, y)
         + tex2D(texRef01, x - 21, y) + tex2D(texRef01, x - 22, y)
         + tex2D(texRef01, x - 23, y) + tex2D(texRef01, x - 24, y)) / 16;
      float average4 = (tex2D(texRef01, x + 9, y) + tex2D(texRef01, x + 10, y)
         + tex2D(texRef01, x + 11, y) + tex2D(texRef01, x + 12, y)
         + tex2D(texRef01, x + 13, y) + tex2D(texRef01, x + 14, y)
         + tex2D(texRef01, x + 15, y) + tex2D(texRef01, x + 16, y)
         + tex2D(texRef01, x + 17, y) + tex2D(texRef01, x + 18, y)
         + tex2D(texRef01, x + 19, y) + tex2D(texRef01, x + 20, y)
         + tex2D(texRef01, x + 21, y) + tex2D(texRef01, x + 22, y)
         + tex2D(texRef01, x + 23, y) + tex2D(texRef01, x + 24, y)) / 16;
      average3 = average3 > average4 ? average3 : average4;
      cfar[tid] = (cfar[tid] >= average3 * threshold) ? cfar[tid] : 0;
   }
}

同时，我使用了6个流，但是计算时间远超传输时间，因此效果很不明显。
其NVP分析结果如下：

我现在有两个想法：
一是，Hyper-Q,可以实现类似这种效果：

但是我不知道怎么做，能不能实现，我的显卡本身是帕斯卡架构，不是开普勒架构。
第二个想法是，参考DSP，使用乒乓操作，像这样：

箭头表示同时进行，这样就不必纠结总的处理时间，只要保证每个步骤的处理时间在200ms以内就行了。
说到这里，其实两种思路是一样的，只不过第一个是并发执行同类内核，第二个是并发执行不同内核，求[名词6]给点意见，不知道怎么做了，我觉得对于P40来说，这种计算量应该不算是规模大的计算，毕竟占用率最高也就9.4%。实在不行，我把代码贴出来，求[名词6]看看哪里可以优化的。

我目前正在尝试使用GPU并行计算处理雷达信号，需要在200ms以内同时完成6批数据的一系列处理。每批的数据量是4096×512，我使用的显卡是Tesla P40.
我的gridDim是（32，16，16），blockDim是（8，8，4）
目前我遇到的问题是：有一个计算环节，我将数据排列成二维的矩阵，需要对每个数据单元，取其上下左右各16个数据，共64个数据求和取平均，计算量很大。为了处理边界问题，我还使用了纹理数组，但即使如此，计算时间还是超出预算....

关键不在于你是否使用了texture。使用了也不一定快.

根据你的代码看，你这等于是native的直接实现。将你这里的tex替换成普通的访问（刚才说了，在Pascal上，用texure不一定总是更快的，实际上，它和普通的数据Cache都一样的，已经合并成了Unified Cache了），看成是普通访存的话，你等于直接的将每个点，反复的取左边，右边，上面，下面的16个点。这种做法并非最优的（从算法的角度）。

该部分的代码如下：
__global__ void CFAR0(float *cfar, float *real)
{
   int bid = blockIdx.x + blockIdx.y*gridDim.x + blockIdx.z*gridDim.x*gridDim.y;
   int tid = bid * blockDim.x * blockDim.y*blockDim.z + threadIdx.x + threadIdx.y * blockDim.x + threadIdx.z*blockDim.x*blockDim.y;
   float x = tid % DistanceDim;
   float y = tid / DistanceDim;
   float average1 = (tex2D(texRef01, x, y - 9) + tex2D(texRef01, x, y - 10)
      + tex2D(texRef01, x, y - 11) + tex2D(texRef01, x, y - 12)
      + tex2D(texRef01, x, y - 13) + tex2D(texRef01, x, y - 14)
      + tex2D(texRef01, x, y - 15) + tex2D(texRef01, x, y - 16)
      + tex2D(texRef01, x, y - 17) + tex2D(texRef01, x, y - 18)
      + tex2D(texRef01, x, y - 19) + tex2D(texRef01, x, y - 20)
      + tex2D(texRef01, x, y - 21) + tex2D(texRef01, x, y - 22)
      + tex2D(texRef01, x, y - 23) + tex2D(texRef01, x, y - 24)) / 16;

   cfar[tid] = (real[tid] >= average1 * threshold) ? real[tid] : 0;
   if (cfar[tid] != 0){
      float average2 = (tex2D(texRef01, x, y + 9) + tex2D(texRef01, x, y + 10)
         + tex2D(texRef01, x, y + 11) + tex2D(texRef01, x, y + 12)
         + tex2D(texRef01, x, y + 13) + tex2D(texRef01, x, y + 14)
         + tex2D(texRef01, x, y + 15) + tex2D(texRef01, x, y + 16)
         + tex2D(texRef01, x, y + 17) + tex2D(texRef01, x, y + 18)
         + tex2D(texRef01, x, y + 19) + tex2D(texRef01, x, y + 20)
         + tex2D(texRef01, x, y + 21) + tex2D(texRef01, x, y + 22)
         + tex2D(texRef01, x, y + 23) + tex2D(texRef01, x, y + 24)) / 16;
      cfar[tid] = (cfar[tid] >= average2 * threshold) ? cfar[tid] : 0;
   }
   if (cfar[tid] != 0){
      float average3 = (tex2D(texRef01, x - 9, y) + tex2D(texRef01, x - 10, y)
         + tex2D(texRef01, x - 11, y) + tex2D(texRef01, x - 12, y)
         + tex2D(texRef01, x - 13, y) + tex2D(texRef01, x - 14, y)
         + tex2D(texRef01, x - 15, y) + tex2D(texRef01, x - 16, y)
         + tex2D(texRef01, x - 17, y) + tex2D(texRef01, x - 18, y)
         + tex2D(texRef01, x - 19, y) + tex2D(texRef01, x - 20, y)
         + tex2D(texRef01, x - 21, y) + tex2D(texRef01, x - 22, y)
         + tex2D(texRef01, x - 23, y) + tex2D(texRef01, x - 24, y)) / 16;
      float average4 = (tex2D(texRef01, x + 9, y) + tex2D(texRef01, x + 10, y)
         + tex2D(texRef01, x + 11, y) + tex2D(texRef01, x + 12, y)
         + tex2D(texRef01, x + 13, y) + tex2D(texRef01, x + 14, y)
         + tex2D(texRef01, x + 15, y) + tex2D(texRef01, x + 16, y)
         + tex2D(texRef01, x + 17, y) + tex2D(texRef01, x + 18, y)
         + tex2D(texRef01, x + 19, y) + tex2D(texRef01, x + 20, y)
         + tex2D(texRef01, x + 21, y) + tex2D(texRef01, x + 22, y)
         + tex2D(texRef01, x + 23, y) + tex2D(texRef01, x + 24, y)) / 16;
      average3 = average3 > average4 ? average3 : average4;
      cfar[tid] = (cfar[tid] >= average3 * threshold) ? cfar[tid] : 0;
   }
}

回到你的代码看，最好的做法是不使用texure/data cache（unified cache），而使用shared memory/warp shuffle。这样可以在不改动算法的情况下，取得更好的效果，容我慢慢说来：

（1）Pascal架构上，L1 cache相比shared memory容量较小，速率较慢，延迟较高。
（2）在此架构的卡上（你的GP102），L1还存在一定情况下失效的问题，此时将直接使用L2。而L2就比Global meory的带宽高一点点。不利于你的性能发挥。
（3）因为一个矩形区域中的两个点（例如点(3,4)和点(3,3)，和点(2,4)）存在大量的重复数据，直接使用shared memory进行手工安排很方便，同时此时，更小的延迟，更大的带宽，最大的容量，一定总是可用，不会失效，往往会带来性能提升。
（4）至于越界，越界不是必须使用texture的，你可以可以手工处理。
（5）当改写成shared memory的时候，你可以进一步的改成使用warp shuffle，将数据存在在寄存器中，因为你的下标规律都是很规则的，可以很方便的改写成warp shuffle。此时可能会进一步的提升性能。

同时，我使用了6个流，但是计算时间远超传输时间，因此效果很不明显。
其NVP分析结果如下：

关于你的是否能用上Hyper-Q(多硬件命令队列）功能，该功能从计算能力5.0+起，已经不再宣传了。因为它已经集成到了所有的卡中（包括家用卡），只要你不使用WDDM驱动（Windows上使用TCC驱动），或者Linux下使用，均会默认启用的。这个并不需要单独处理。

箭头表示同时进行，这样就不必纠结总的处理时间，只要保证每个步骤的处理时间在200ms以内就行了。
说到这里，其实两种思路是一样的，只不过第一个是并发执行同类内核，第二个是并发执行不同内核，求[名词6]给点意见，不知道怎么做了，我觉得对于P40来说，这种计算量应该不算是规模大的计算，毕竟占用率最高也就9.4%。实在不行，我把代码贴出来，求[名词6]看看哪里可以优化的。

将代码拆分成多个处理部分的kernel（例如拆分成处理左边，处理右边，处理上边，处理下面）并无意义，甚至因为中间状态需要保存到global memory, 以便在多个kernel间反复的写入+读取，反而可能会降低性能。

你的根本问题在于，你没有选好算法。如果你要坚持目前的实现，则“改善性”的方案就是刚才说的（shared memory和shuffle），但是存在如下数学/算法上的变通：

考虑到点(x=100,y=100)和点(x=101, y=100), 他们都计算自己左侧的16个点的时候，实际上，其中15个点的加法的值是一样的，你这样做，等于冗余计算了很多倍。类似的还有上面，下面，右面。

这问题实际上是一个moving avarage的问题，如果点p0 = (x,y)自己左侧的16个点的平均值是a0的话，则(x+1,y)左侧16个点的平均值可以通过a0 + (val(x,y) - (val(x-16,y))得到，这样每次只要只要干掉一个旧点，更新上一个心点，即可得到值，这样对访存的要求就大量降低了（无论你是global memory直接，还是texture，还是shared memory），对计算量的要求也降低了很多了。这样才算真正的优化。

这的确是一个滑动的窗口，无论是上下还是左右，滑动的时候，每次出去一个旧值，就进来一个新值。

我建议你按照这个新算法进行改动，但具体映射到CUDA上的实现是你的问题（例如，如何有效的利用旧的a0值？如何有效的纵向（看起来）访存）。但这个毕竟比你的老算法有相当大的计算量和访存上上的优势，我不建议你坚持老算法（如果要坚持，请看之前的建议）。

中午打了很大的一段话，上千字吧，可惜没有了。这次打一段发一段吧。

非常感谢您给了我一个明确的方向，其实我有往这方面想，但是没有专注研究过，总觉得要分类、要判断就会变麻烦、变低效，往深了说就是怕麻烦 ::)我使用纹理也不是因为觉得它快，而是因为它不用讨论边界问题，所有单元的动作都一致，没别的意思。