前言

《并行编程》系列是学习《Intro to Parallel Programming》过程中所做的笔记记录以及个人一些所思所想。

并行通信

并行计算需要解决的一个问题就是，如何解决线程之间的协同工作（Working together）问题。而协同工作的关键则是通信（Communication）。

CUDA 的通信发生在内存当中，例如，多个线程可能需要从同一个内存地址中读取数据；也可能出现多个线程需要同时向同一个地址写入数据；可能出现多个线程之间需要交换数据。

并行通信模式（Parallel Communication Patterns）

并行通信存在多种模式，通信模式反映了线程所执行的任务与内存之间的映射关系。这里将分别介绍五种通信模式：

1. Map
2. Gather
3. Scatter
4. Stencil
5. Transpose

Map

Map: Tasks read from and write to specific data elements.

Map 模式下，每个线程将从内存的特定地址中读取数据进行处理，然后再写入特定的地址中，它的输入与输出具有严格的一对一的关系。

Map 模式在 GPU 中非常高效，在 CUDA 中也能很容易通过有效的方式表达。

但是 Map 比较不灵活，能处理的问题有限。

Gather

现在假设需要求取3个数据的平均值，那么在 Gather 模式下，每个线程将从内存中的三个位置读取数据，然后将这三个数取平均，写入指定的内存中。

这一模式可用于涉及到滤波器的一系列操作。

Scatter

Scatter: Tasks compute where to write output.

与 Gather 模式下，多个输入一个输出相反，Scatter 模式是一个输入多个输出。

另外在同时写入多个输出的时候将出现冲突的问题，这将在后续进行讨论。

Stencil

Stencil: Tasks read input from a fixed neighborhood in an array.

常用的模板有：

1. 2D von Neumann

   
   

2. 2D Moore

   
   

3. 3D von Neumann

   
   

看到这里，可能有人会对 Stencil 和 Gather 产生疑惑。咋看之下，两者确实非常相似，但是 Stencil 模式中，要求每个线程都严格执行相同的模板，但是 Gather 模式却没有这个限制，因此，比如说，在 Gather 模式中就可以按线程索引的奇偶不同，给线程分配不同的操作任务。

Transpose

Transpose: Tasks re-order data elements in memory.

对于一张图像，其数据在内存中的存储的方式如下：

但是在某些情况下，可能需要将图像转置。

通常在涉及到数组运算、矩阵运算和图像操作的时候会需要使用到 Transpose，但是 Transpose 也适用于其它数据结构。

比如定义了一个结构体 foo ，然后创建一个该结构的结构数组（AoS），如果想将该结构数组变换成数组结构（SoA），也可以通过 Transpose 实现。

总结

上图总结了并行计算的七种计算模式，除了之前介绍的五种模式以外，还有两种更加基础的模式将在接下来进行介绍。

GPU

程序员眼中的 GPU

程序员在并行编程中所要做的就是，创建内核（C/C++函数）用来处理具体的任务。内核由许多线程（完整执行一段处理程序的通路）组成，图中的线程都采用曲线绘制，其原因是，每个线程的具体通路可能不相同（即每个线程所执行的运算不相同）。

多个线程将组成线程块，一个线程块内的多个线程负责协同处理一项任务或者子任务。

上图中，程序首先启动了一个内核 foo，等到其中所有的线程都运行完了之后，结束内核。然后又启动了内核 bar，可以注意到，一个内核中所具有的线程块，以及每个线程块中的线程数是可以自己配置的参数。

线程块与 GPU

GPU 中包含有许多的流处理器（Streaming Multiprocessor, SM），不同的 GPU 包含有不同数量的流处理器，并且流处理器数量也是衡量 GPU 性能的一项重要指标。

一个流处理器中包含有多个简单的处理器和内存。

当你的程序创建了内核之后，GPU 将为内核中的线程块分配流处理器，每个线程块被分配给一个流处理器，然后这些流处理器以并行的方式进行运行。

注意：一个流处理器上允许运行多个线程块，但是一个线程块只允许被分配给一个流处理器运行。

CUDA 特征

CUDA 不具备的特征

CUDA 对于内核中的线程块要何时运行、该如何分配流处理以及有多少线程块需要同时运行等细节没有进行任何的控制，这些分配问题都交给 GPU 进行控制。这么做的好处有：

• 硬件将可以更加高效地执行计算
• 当一个线程块执行完成之后，当前的流处理器马上又可以任意执行下一个线程块
• 更高的扩展性。因为流处理器的分配交由硬件控制，所以程序可以很好地在具有不同流处理器数量的设备上进行移植。

但是 CUDA 的这种做法也将导致一些后果：

• 对于某一线程块将在哪个流处理器上运行无法做出任何预测
• 线程块之间没有通信交流。如果线程块 x 的输入依赖于线程块 y 的输出，而 y 已经完成执行并且退出，这将导致 x 的计算出现问题。这种现象称为“dead lock”
• 线程块中的线程不能永远执行（比如，死循环），因为它需要在执行完成之后释放流处理器资源，以便于其它线程块可以使用

CUDA 具备的特征

CUDA 在程序运行的时候，能够保证两点：

• 同一个线程块上的所有线程将同时在同一个流处理器中运行。
• 下一个内核中的线程块必须等待当前内核中的所有块运行完成之后，才能运行。
  • 比如说，程序依次定义了两个内核 foo 和 bar ，bar 中的线程块只有等到 foo 中的所有线程块都运行完之后才能开始运行。

GPU 内存模型

每个线程都拥有一个局部内存（Local memory），这就好像局部变量一样，只有对应的线程才能访问。

然后，线程块也有一块对应的共享内存（Shared memory）。共享内存只能被对应线程块内的线程进行访问。

另外还有具有全局内存（Global memory）。不仅内核中的所有线程可以访问它，不同内核也可以进行访问。

前边介绍的局部内存、共享内存和全局内存都是属于 GPU 内部的内存。上图展示了，CPU 的线程启动了 GPU ，然后将主机内存（Host memory）中的数据拷贝到 GPU 的全局内存中，以便于 GPU 内核线程可以访问这些数据。另外 GPU 内核线程也可以直接访问主机内存，这一点将在后边介绍。

同步

通过共享内存和全局内存，线程之间可以互相访问彼此的计算结果，这也意味着线程间可以进行协同计算。但是这样也存在着风险， 如果一个线程在另一个线程写入数据之前就读取了数据怎么办？因此线程之间需要同步的机制，来避免这种情形出现。

事实上，同步问题是并行计算的一个最基本的问题。而解决同步问题的一个最简单方法则是屏障（Barrier）。

Barrier: Point in the program where threads stop and wait. When all threads have reached the barrier, they can proceed.

屏障语句是 __syncthreads() 。

编程模型

现在，可以重新构建一下编程模型。我们拥有线程和线程块，并且在线程块内，可以创建屏障用于同步线程。事实上，如果一个程序中创建了多个内核，内核之间默认具有隐性的屏障，这使得不会出现多个内核同时运行的情况。

然后再将之前介绍的内存模型添加进来，便得到了 CUDA 。

因此，CUDA 的核心就是层级计算结构。从线程到线程块再到内核，对应着内存空间中的局部内存、共享内存和全局内存。

编写高效的 CUDA 程序

这里将首先从顶层的策略上介绍如何编写高效的 CUDA 程序。

首先需要知道的是 GPU 具有非常惊人的计算能力，一个高端的 GPU 可以实现每秒超过 3 万亿次的数学运算（3 TFLOPS/s）。但是如果一个 CUDA 程序的大多数时间都花费在了等待内存的读取或写入操作的话，这就相当浪费计算能力。所以要编写高效的 CUDA 程序的第一点是——最大化计算强度。

计算强度表达为每个线程计算操作时间除以每个线程在的访存时间。所以要最大化计算强度，就可以通过最大化分子和最小化分母来实现。然而由于计算操作时间主要受具体算法的计算量限制，所以为了最大化计算强度主要从最小化访存时间入手。

最小化访存时间

要最小化访存时间的一种方式就是，将访问频率更高的数据移动到访问速度更快的内存中。

在之前的介绍当中已经了解了 GPU 线程可以访问四种类型的内存，其中最快就是局部内存。

局部内存

局部变量的定义是最简单的。

对于上图的内核代码，变量 f 与参数 in 都将存储于局部内存中。

共享内存

要定义存储于共享内存中的变量，需要在变量定义语句前加一个 __shared__ 关键字进行修饰。定义于共享内存中的变量可以被同一个线程块中的所有线程所访问，其生存时间为线程块的生存时间。

全局内存

全局的内存访问要稍微麻烦些，但是可以通过指针的机制来实现。

这里传入内核的参数被定义成一个指针，而这个指针恰恰指向的是全局内存区域。

然后在 CPU 的代码部分，首先创建了一个长度为 128 的浮点数数组 h_arr，它将存储于主机内存中（这里通过前缀 h_ 表明当前变量运行于 HOST 中），然后定义了一个指向 GPU 全局内存的指针 d_arr，并通过 cudaMalloc 函数为 d_arr 分配全局存储区域。

最小化访存时间的另一个方法是使用合并全局内存访问（Coalesce global memory accesses）。

单一线程在访问内存时具有一个特性，就是即使该线程只需要使用到内存中的一小部分，但是程序也会从内存中读取一段连续的内存块。因此，如果此时恰好有其它线程也在使用该内存块中的数据，内存块就得到复用，从而节省再次读取内存的时间。

所以如果多个线程同时读取或者写入连续的全局内存位置，此时 GPU 的效率的是最高的，而这种访问模式被称为合并（Coalesced）。

但是当多个线程所访问的全局内存位置不连续或者甚至随机的时候，此时 GPU 便无法继续保持高效，因为很可能需要分别读取全局内存中的多个块，这样就增加了访存时间。

相关性问题（Related problem）

Related problem: lots of threads reading and writing same memory locations

当多个线程同时参与到对同一块内存地址的读写操作时，将引发冲突从而导致错误的计算结果，这便是相关性问题。

解决该相关性问题的一个方法是使用原子内存操作（Atomic memory operations）。

原子内存操作

CUDA 提供了若干个原子内存操作函数，通过这些函数可以以原子操作的方式访问内存，也就是某一时刻内存中的特定地址只能被单一线程所读写，从而避免了相关性问题。

常见的原子内存操作：

• atomicAdd()，原子相加
• atomicMin()，原子最小值
• atomicXOR()，原子异或
• atomicCAS()，比较并且交换（Compare-and-Swap）

说明：这些原子内存操作函数的实现借助了硬件来实现原子操作，这里将不进行介绍。

但是这些原子操作也存在一些局限性。

• 只支持某些特定的操作（比如，支持加、减、最小值和异或等，不支持求余、求幂等操作）和数据类型（主要支持整数）。
• 没有顺序限制。尽管使用了原子操作，但是关于线程执行顺序的问题依然没有定义。
  • 由于浮点数精度问题，这将导致浮点数运算出现非关联现象（Non-associative）。具体来说就是可能出现 (a + b) + c != a + (b + c)，比如，当 a = 1, b = 10^99, c= 10^-99 时。
• 串行化线程内存访问。原子操作的实现并没有使用什么神奇的魔法，它仅仅只是串行化了线程对同一个内存地址的访问，所以这将减慢整体的计算速度。

线程发散

前边已经介绍过了，要使得 CUDA 程序高效的一个关键点是——最大化计算强度。然后另外一个关键点是——避免线程发散（Thread divergence）。

线程发散指的是，比如说当内核代码中出现条件语句时，线程运行到条件语句处，可能有些线程符合条件，而有些线程不符合条件，此时它们就会发散开，形成两条路径，然后在条件语句块结束之后再次聚合到同一条路径上。

不仅仅只有条件语句才会导致线程发散，循环语句也可能导致。

举个不太恰当的例子，在这个内核代码中有一个循环，循环的次数是当前线程的索引。

所以线程的执行路径如上图，如果以时间为横轴绘制线程运行图则如下图。

由于硬件倾向于同时执行完线程，所以当线程索引小的线程完成循环之后，它还会继续等待其它线程完成循环，直至所有线程都完成循环之后，这些线程才会继续执行循环块之后的代码。因此，这里除了最后一个线程充分利用了时间进行运算以外，其它线程均无法有效利用时间。而这也就是为什么要避免线程发散的原因。

总结

本节内容小结：

• 通信模式
  • gather, scatter, stencil, transpose
• GPU 硬件与编程模型
  • 流处理器，线程，线程块
  • 线程同步
  • 内存模型（局部，共享，全局，主机），原子操作
• 高效 GPU 编程
  • 减少访存花销（使用更快的内存，合并全局内存访问）
  • 避免线程发散

课堂作业

本次的课堂作业是实现图像模糊，思路相对较简单。唯一需要注意的是边界情况的取值。因为当 filter 的中心位于图像边界的时候，它的周围像素会出现超出图像的现象，这里需要进行判断。

课程作业完成代码：
https://github.com/un-knight/cs344-parallel-programming