CUDA C编程（十）核函数可达到的带宽

  在分析核函数性能时，需要注意内存延迟，即完成一次独立内存请求的时间；内存带宽，即SM访问设备内存的速度，它以每单位时间内的字节数进行测量。我们通常尝试使用以下两种方法来改进核函数的性能：1.通过最大化并行执行线程束的数量来隐藏内存延迟，通过维持更多正在执行的内存访问来达到更好的总线利用率；2.通过适当的对齐和合并内存访问来最大化内存带宽效率。然而往往当前问题的本质就是有一个不好的访问模式，接下来我们将利用一个矩阵转置的例子学习如何使用各种优化手段来调整核函数的带宽。

内 存 带 宽
  大多数核函数对内存带宽非常敏感，也就是说它们有内存带宽的限制。因此，在调整核函数时需要注意内存带宽的指标。全局内存中数据的安排方式，以及线程束访问该数据的方式对带宽有显著影响。一般有两种类型的带宽：理论带宽和有效带宽。理论带宽时当前硬件可以实现的绝对最大带宽。有效带宽是核函数实际达到的带宽，它是测量带宽，可以用下列公式计算：

矩 阵 转 置 问 题
  矩阵转置是线性代数中一个基本问题，虽然是基本问题，但却在许多应用中被使用。矩阵的转置意味着每一列与相应的一行进行互换，如下图所示：

  以下是基于主机实现的使用单精度浮点值的错位转置算法。假设矩阵存储在一个一维数组中。通过改变数组索引值来交换行和列的坐标，可以很容易得到转置矩阵。

void transposeHost(float *out, float *in, const int nx, const int ny)
{
   for(int iy = 0; iy < ny; ++iy)
   {
      for(int ix = 0; ix < nx; ++ix)
      {
         out[ix * ny + iy] = in[iy * nx + ix];
      }
   }
}

  在这个函数中有两个用一维数组存储的矩阵：输入矩阵in和转置矩阵out。矩阵维度被定义为nx行ny列。可以用一个一维执行转置操作，如下图所示：

  观察输入和输出布局，我们会注意到：读：通过原矩阵的行进行访问，结果为合并访问；写：通过转置矩阵的列进行访问，结果为交叉访问。其中交叉访问是使GPU性能变得最差的内存访问模型，但是，在矩阵转置操作中这是不可避免地。接下来开始侧重于使用两种转置核函数来提高带宽的利用率：一种是按行读取按列存储，另一种是按列读取按行存储。下图是两种方法的实现，如果禁用一级缓存加载，那么这两种实现的性能在理论上是相同的。但是如果启用一级缓存，那么第二种实现的性能表现得更好。按列读取操作是不合并的(因此带宽会被浪费在未被请求的字节上)，将这些额外的字节存入一级缓存意味着下一个读操作可能会在缓存上执行而不再全局内存上执行。因为写操作不在一级缓存中缓存，所以对按列执行写操作的例子而言，任何缓存都没有意义。

为转置核函数设置性能的上限和下限
  在执行矩阵转置核函数之前，可以先创建两个拷贝核函数来粗略计算所有转置核函数性能的上限和下限：通过加载和存储行来拷贝矩阵(上限)，这样将模拟执行相同数量的内存操作作为转置，但是只能使用合并访问；通过加载和存储列来拷贝矩阵(下限)，这样将模拟相同数量的内存操作作为转置，但是只能使用交叉访问。核函数的实现如下：

__global__ void copyRow(float *out, float *in, const int nx, const int ny)
{
   unsigned int ix = blockDim.x * blockIdx.x + threadIdx.x;
   unsigned int iy = blockDim.y * blockIdx.y + threadIdx.y;
   if(ix < nx && iy < ny)
   {
       out[iy * nx + ix] = in[iy * nx + ix];
   }
}

__global__ void copyCol(float *out, float *in, const int nx, const int ny)
{
   unsigned int ix = blockDim.x * blockIdx.x + threadIdx.x;
   unsigned int iy = blockDim.y * blockIdx.y + threadIdx.y;
   if(ix < nx && iy < ny)
   {
       out[ix * ny + iy] = in[ix * ny + iy];
   }
}

朴素转换：读取行和列
  基于行的朴素转置核函数是基于主机实现的，这种转置按行加载按列存储：

__global__ void transposeNativeRow(float *out, float *in, const int nx, const int ny)
{
   unsigned int ix = blockDim.x * blockIdx.x + threadIdx.x;
   unsigned int iy = blockDim.y * blockIdx.y + threadIdx.y;

   if(ix < nx && iy < ny)
   {
      out[ix * ny + iy] = in[iy * nx + ix];
   }
}

  通过互换读索引和写索引，就生成了基于列的朴素转置核函数，这种转置按列加载按行存储:

__global__ void transposeNativeCol(float *out, float *in, const int nx, const int ny)
{
   unsigned int ix = blockDim.x * blockIdx.x + threadIdx.x;
   unsigned int iy = blockDim.y * blockIdx.y + threadIdx.y;
   if(ix < nx && iy < ny)
   {
       out[iy * nx + ix] = in[ix * ny + iy];
   }
}

  结果表明，使用NativeCol方法比NativeRow方法性能表现得更好。导致这种性能提升的原因时在缓存中执行了交叉读取。即使通过某一方式读取一级缓存中的数据没有都被这次访问使用到，这些数据仍留在缓存中，在以后的访问过程中可能发生缓存命中。
展开转置：读取行和读取列
  接下来将利用展开技术来提高转置内存带宽的利用率，其中，展开的目的是为每个线程分配更独立的任务，从而最大化当前内存需求。以下是一个展开因子为4的基于行的实现。这里引入了人两个新的数组索引：一个用于行访问，另一个用于列访问。

__global__ void transposeUnroll4Row(float *out, float *in, const int nx, const int ny)
{
   unsigned int ix = blockDim.x * blockIdx.x * 4 + threadIdx.x;
   unsigned int iy = blockDim.y * blockIdx.y + threadIdx.y;
   unsigned int ti = iy * nx + ix;
   unsigned int to = ix * ny + iy;
   if(ix + 3 * blockDim.x < nx && iy < ny)
   {
      out[to] = in[ti];
      out[to + ny * blockDim.x] = in[ti + blockDim.x];
      out[to + ny * 2 * blockDim.x] = in[ti + 2 * blockDim.x];
      out[to + ny * 3 * blockDim.x] = in[ti + 3 * blockDim.x];
   }
}

  使用相似的展开交换读索引和写索引产生一个基于列的实现：

__global__ void transposeUnroll4Row(float *out, float *in, const int nx, const int ny)
{
   unsigned int ix = blockDim.x * blockIdx.x * 4 + threadIdx.x;
   unsigned int iy = blockDim.y * blockIdx.y + threadIdx.y;
   unsigned int ti = iy * nx + ix;
   unsigned int to = ix * ny + iy;
   if(ix + 3 * blockDim.x < nx && iy < ny)
   {
      out[ti] = in[to];
      out[ti + 1 * blockDim.x] = in[to + ny * 1 * blockDim.x];
      out[ti + 2 * blockDim.x] = in[to + ny * 2 * blockDim.x];
      out[ti + 3 * blockDim.x] = in[to + ny * 3 * blockDim.x];
   }
}

对角转置：读取行与读取列
  当启用一个线程块的网格时，线程块会被分配给SM。编程模型抽象可能用一个一维或二维布局来表示该网格，但是从硬件的角度来看，所有块都是一维的。每个线程块都有其唯一标识符bid，它可以用网格中的线程块按行顺序优先顺序计算得出：int bid = blockIdx.y * gridDim.x + blockIdx.x;下图所示是一个4×4的线程块网格，它包含了每个线程块的ID。当启用一个核函数时，线程块被分配给SM的顺序由块ID来确定。一旦所有的SM被完全占用，所有剩余的线程块都保持不变直到当前的执行被完成。一旦一个线程块执行结束，将为该SM分配另一个线程块。由于线程块完成的速度和顺序是不确定的，随着内核进程的执行，起初通过bid相连的活跃线程块会变得不太连续了。

  尽管无法直接调控线程块的顺序，但我们可以灵活的使用块坐标blockIdx.x和blockIdx.y。上图说明了笛卡尔坐标系下的块坐标。下图展示了一个标识blockIdx.x和blockIdx.y的不同方法：使用对角块坐标系。

  对角坐标系用于确定一维线程块的ID，但对于数据访问，仍需要使用笛卡尔坐标系。因此，当用对角坐标表示块ID时，需要将对角坐标映射到笛卡尔坐标中，以便可以访问到正确的数据块。对于一个方阵来说，这个映射可以通过以下方程式计算得出：block_x = (blockIdx.x + blockIdx.y) % gridDim.x; block_y = blockIdx.x;这里的blockIdx.x和blockIdx,y为对角坐标。block_x和block_y是它们对应的笛卡尔坐标，基于行的矩阵转置核函数使用如下所示的对角坐标。在核函数的起始部分包含了从对角坐标到笛卡尔坐标的映射计算，然后使用映射的笛卡尔坐标(block_x,block_y)来计算线程索引（ix,iy），这个对角转置核函数会影响线程块分配数据块的方式。下面的核函数使用了对角线程块坐标，它借助合并读取和交叉写入实现了矩阵的转置：

__global__ void transposeDiagonalRow(float *out, float *in, const int nx, const int ny)
{
   unsigned int blk_y = blockIdx.x;
   unsigned int blk_x = (blockIdx.x + blockIdx.y) % gridDim.x;
   unsigned int ix = blockDim.x * blk_x + threadIdx.x;
   unsigned int iy = blockDim.y * blk_y + threadIdx.y;

   if(ix < nx && iy < ny)
   {
      out[ix * ny + iy] = in[iy * nx + ix];
   }
}

  使用基于列的对角坐标的核函数如下所示：

__global__ void transposeDiagonalRow(float *out, float *in, const int nx, const int ny)
{
   unsigned int blk_y = blockIdx.x;
   unsigned int blk_x = (blockIdx.x + blockIdx.y) % gridDim.x;
   unsigned int ix = blockDim.x * blk_x + threadIdx.x;
   unsigned int iy = blockDim.y * blk_y + threadIdx.y;

   if(ix < nx && iy < ny)
   {
      out[iy * nx + ix] = in[ix * ny + iy];
   }
}

  通过使用对角坐标系来修改线程块的执行顺序，这使得基于行的核函数1性能得到了大大提升。但是，基于列的核函数在使用笛卡尔块坐标系仍比使用对角块坐标系表现得更好。对角块核函数的实现可以通过展开块得到更大的提升，但是这种实现不像使用基于笛卡尔坐标系的核函数那样简单直接。这种性能的提升的原因与DRAM的并行访问有关，发送给全局内存的请求由DRAM分区完成。设备内存中的连续的256字节区域被分配到连续的分区。当使用笛卡尔坐标将线程块映射到数据块时，全局内存访问可能无法均匀地被分配到整个DRAM分区中，这时就可能发生“分区冲突”的现象。发生分区冲突时，内存请求在某些分区中排队等候，而另一些分区一直未被调用。因为对角坐标映射造成了从线程块到待处理数据块的非线性映射，所以交叉访问不太可能会落入到一个独立的分区中，并且会带来性能的提升。
使用瘦块来增加并行性
  增加并行性的方式是调整块的大小，块大小的测试结果如下表所示：

  目前最佳的块大小为（8，32），尽管它与大小为(16,16)的块显示了相同的并行性，但这种性能的提升是由"瘦的"块(8,32)带来的，如下图所示：

  通过增加存储在线程块中连续元素的数量，“瘦”的块可以提高存储操作的效率。