23、并行算法：原理、模型与问题求解

并行算法：原理、模型与问题求解

1. 计算机系统架构分类

如今，耗时的计算任务以及多处理器计算机系统的动态发展，使得开发和分析并行算法的问题成为焦点。由于存在各种各样的计算机系统架构，因此有必要根据各种参数对它们进行分类。历史上，Flynn 提出了一种基于多指令流和数据流标准划分架构的方法。

流被定义为处理器执行或处理的指令或数据序列。从这个角度来看，程序为处理器提供要执行的指令流，用于处理的数据也以流的形式出现。根据 Flynn 的分类法，假设指令流和数据流是独立的，因此计算系统可分为以下几类：
1. 单指令流，单数据流（SISD） ：具有单核处理器的标准计算机具有这些特性，它们一次只能执行一个操作。
2. 单指令流，多数据流（SIMD） ：在这类系统中，相同的操作会同时在不同的数据上执行。例如，向量计算机系统就属于这一类，其中一条指令可以对一组数据元素执行操作。
3. 多指令流，单数据流（MISD） ：尽管这种方法缺乏实际意义，但 MISD 机器对于某些高度专业化的任务可能有用。
4. 多指令流，多数据流（MIMD） ：MIMD 系统是 Flynn 分类中最广泛和多样的类别。大多数现代多处理器计算机系统都属于这一类。

这些计算机系统类别的特性总结在下表中：
| 计算机系统 | 单数据流 | 多数据流 |
| — | — | — |
| 单指令流 | SISD
a1 + b1 | SIMD
a1 + b1
a2 + b2
a3 + b3 |
| 多指令流 | MISD
a1 + b1
a1 - b1
a1 * b1 | MIMD
a1 + b1
a2 - b2
a3 * b3 |

此外，对 Flynn 的分类法有一个广泛使用的细化，根据计算机系统中内存组织的方式，MIMD 类别又可分为多处理器（统一内存访问，UMA）和多计算机（无远程内存访问，NORMA）。多计算机中处理器之间的交互通过消息传递机制执行。

2. PRAM 模型

在之前的研究中，我们使用了随机访问机（RAM）模型。RAM 系统由处理器、内存访问设备（系统总线）和由有限数量的单元组成的内存构成。处理器依次执行程序的指令，执行主要的算术和逻辑运算，并在内存中读写数据。假设每条指令在固定的时间间隔内执行。处理器的随机操作包括三个阶段：
1. 从内存读取数据到其寄存器 ri 之一（1 ≤ i ≤ N）。
2. 对其寄存器的内容执行算术或逻辑运算。
3. 将数据从寄存器 rj（1 ≤ j ≤ N）写入某个内存单元。

假设上述三个步骤的执行时间为 Θ(1)。

并行随机访问机（PRAM）是最广泛使用的并行计算机系统模型之一。PRAM 结合了 p 个处理器、共享内存和一个将程序指令传输到处理器的控制设备。PRAM 的一个重要特征是系统中任何处理器对随机内存单元的访问时间有限。与 RAM 一样，算法步骤对应于三个处理器操作：
1. 处理器 Pi 从第 j 个内存单元读取数据。
2. 处理器 Pi 对其寄存器的内容执行算术或逻辑运算。
3. 将数据写入第 k 个内存单元。

同样，一个算法步骤的执行时间为 Θ(1)。

当两个或多个处理器同时访问同一个内存单元时，会导致访问冲突，可分为读冲突和写冲突。如果多个处理器试图从一个单元读取数据，有两种操作选项：
- 独占读（ER） ：在给定时间内只允许一个处理器读取，否则程序出错。
- 并发读（CR） ：访问同一内存单元的处理器数量不受限制。

如果多个处理器试图在一个地址写入数据，也有两种选项：
- 独占写（EW） ：在特定时刻只允许一个处理器写入给定单元。
- 并发写（CW） ：多个处理器可以同时访问单个内存单元。

在并发写的情况下，处理器遵循的写入规则有以下几种：
- 通用值记录 ：准备写入单个内存单元的所有处理器必须记录一个通用值，否则写入指令被视为错误。
- 随机选择 ：随机选择执行写入的处理器。
- 优先写入 ：为每个竞争处理器分配一定的优先级，如计算值，只保留来自预先确定优先级（如最低）的处理器的值。
- 混合选择 ：所有处理器提供写入的值，通过某种操作创建结果（如值的总和、最大值等），然后进行记录。

PRAM 根据冲突解决方法的分类如下：

graph LR
    classDef process fill:#E5F6FF,stroke:#73A6FF,stroke-width:2px
    A(PRAM):::process --> B(EREW):::process
    A --> C(CREW):::process
    A --> D(ERCW):::process
    A --> E(CRCW):::process

EREW 系统在处理内存单元时存在显著限制。另一方面，由于同时访问某个内存段的核心数量有限，从技术角度来看，具有大量处理器的 CREW、ERCW、CRCW 系统难以构建。不过，有一个重要且有些意外的结果，即可以在按照 EREW 原则构建的系统上模拟 CRCW 机器的工作。

仿真定理 ：设 CRCW 机器的算法使用 p 个处理器在时间 T(N) 内解决一个规模为 N 的特定问题。那么，存在一个用于 EREW 系统的相同问题的算法，使用 p 个处理器可以在时间 O(T(N) log₂N) 内执行（PRAM 的内存容量必须增加 O(p) 倍）。

与 RAM 模型不同，多处理器计算机系统算法复杂度的主要衡量标准是算法的执行时间。我们引入以下符号：T₁(N) 是顺序算法解决复杂度由参数 N 估计的问题所需的时间；Tp(N) 是并行算法在具有 p 个处理器（p > 1）的机器上所需的时间。由于从 RAM 的定义可知，每个操作都需要一定的时间，因此 T₁(N) 与所使用算法中的计算操作数量成正比。

注意，当 p → ∞ 时，算法的执行时间达到最小值。具有无限数量可用处理器的假设计算机系统称为超计算机。超计算机算法的渐近复杂度表示为 T∞(N)。

为了分析并行算法，广泛使用了加速比、效率和成本的概念。
- 加速比 Sp(N) ：使用具有 p 个处理器的并行算法获得的加速比 Sp(N) 等于 Sp(N) = T₁(N) / Tp(N)。这是与最佳顺序算法相比生产率提高的度量。加速比越大，单处理器系统和多处理器系统上算法运行时间之间的问题解决时间差异就越大。
- 效率 Ep(N) ：特定并行算法使用处理器的效率 Ep(N) 等于 Ep(N) = T₁(N) / (pTp(N)) = Sp(N) / p。
- 成本 Cp(N) ：成本 Cp(N) 定义为并行解决问题的时间与使用的处理器数量的乘积，即 Cp(N) = pTp(N)。成本最优算法的特点是成本与最佳顺序算法的复杂度成正比，此时 Cp(N) / T₁(N) = Θ(1)。

加速比和效率的极限值分别为 Sp = p 和 Ep = 1。当能够在所有处理器上均匀分配计算，并且在程序运行阶段不需要额外的操作来提供处理器之间的通信和组合结果时，可实现 Sp 的最大值。试图通过改变处理器数量来增加加速比会降低 Ep 的值，反之亦然。使用单个处理器（p = 1）时可实现最大效率。

对于 PRAM 模型中的任何算法，都有可能将其修改为适用于处理器数量较少的系统，这就是 Brent 引理：设并行算法 A 为解决某个问题在 RAM 上执行的时间为 T₁，在超计算机上执行的时间为 T∞。那么，存在一个用于解决该问题的算法 A’，在具有 p 个处理器的 PRAM 上执行的时间为 T(p)，其中 T(p) ≤ T∞ + (T₁ - T∞) / p。

3. 求和与部分和问题

考虑计算数组元素 list[1..N]（其中 N = 2^m，m 为正整数）之和的问题。可以通过构建一个完全二叉树来获得解决此问题的并行算法。树的叶子对应于数组元素 list[1..N]，内部节点放置其两个上层节点值的总和，最终在最低层得到所需的总和。这种求和序列的一个特点是树的第 i 层的加法操作与该层的其他操作相互独立。

例如，当 N = 8 时，求和计算如下：

list[1]  list[2]  list[3]  list[4]  list[5]  list[6]  list[7]  list[8]
  [1..2]    [3..4]    [5..6]    [7..8]
    [1..4]        [5..8]
        [1..8]

在本章的程序文本中，使用了保留字 parallel。“parallel for i := 1 to N do …” 表示关于变量 i 的循环迭代应分布到可用的计算机节点上。第一个处理器 P₁ 使用 i = 1 执行循环体的指令，第二个处理器 P₂ 使用 i = 2，依此类推。除非另有说明，假设存在足够数量的处理器 p ≥ N 以实现完全的循环并行化。如果 p < N，算法仍然正确，因为下一个值 i = p + 1 再次由第一个处理器处理，i = p + 2 由第二个处理器处理，以此类推，这种计算负载分配过程称为块循环分配。

以下是计算长度为 N = 2^m 的数组元素并行求和的算法：

const N = 8;
type Mass = array[1..N] of integer;
function parSum(list: Mass): integer;
    // 并行求和算法
    var i, j: integer;
begin
    for i := 1 to Log2(N) do
        parallel for j := N downto 2 step 2^i do
            list[j] := list[j - 2^(i - 1)] + list[j];
    result := list[N];
end;

对于顺序版本，需要 N/2 + N/4 + … + 1 = N - 1 次加法操作。对于具有 p = N/2 个处理器的多处理器计算机系统，找到总和的算法复杂度与“操作 - 操作数”有向图的高度成正比，即 Tp(N) = Θ(log₂N)。因此，加速比 S(N) = Θ(N / log₂N)，效率 E(N) = Θ(log₂⁻¹N)。

部分和计算问题是根据指定的值 list[i]（1 ≤ i ≤ N）确定序列 Sj = ∑(i = 1 to j) list[i] 对于所有 1 ≤ j ≤ N 的值。许多实际问题都可以归结为这个问题，例如确定位串中连续 1 的长度以及连续背包问题。当 j = N 时，就得到了前面分析的数组元素求和问题。

标准顺序算法需要执行 N - 1 次加法操作。使用具有 p = N 个处理器的 PRAM，可以在时间 Θ(log₂N) 内得到结果，通过对数组及其右移 2^(i - 1) 位置的副本进行求和操作，如 N = 8 时的示例：

list[1]  list[2]  list[3]  list[4]  list[5]  list[6]  list[7]  list[8]
  [1]      [1..2]    [2..3]    [3..4]    [4..5]    [5..6]    [6..7]    [7..8]
  [1]      [1..2]    [1..3]    [1..4]    [2..5]    [3..6]    [4..7]    [5..8]
  [1]      [1..2]    [1..3]    [1..4]    [1..5]    [1..6]    [1..7]    [1..8]

以下是计算部分和的算法：

const N = 8;
type Mass = array[1..N] of integer;
procedure parPrefix(list: Mass; var prefixList: Mass);
    // 计算数组元素的部分和
    var i, j: integer;
    temp: Mass;
begin
    for i := 1 to N do prefixList[i] := list[i];
    for i := 1 to Log2(N) do
    begin
        parallel for j := 1 to N do
            temp[j] := prefixList[j];
        parallel for j := N downto 2^(i - 1) + 1 do
            prefixList[j] := temp[j - 2^(i - 1)] + temp[j];
    end;
end;

由于算法第九行中对变量 i 的每个循环迭代中，所有 Θ(N) 次加法操作都是并行执行的，与标准顺序算法相比，加速比 S(N) = Θ(N / log₂N)，效率 E(N) = Θ(log₂⁻¹N)。

将第十行替换为 “for i := 1 to Ceil(Log2(N)) do” 可以使算法 parPrefix 处理任意大小的数组 N = 1, 2, 3, …。算法 parPrefix 的成本是 O(N log₂N)，超过了基本顺序算法的成本。如果将处理器数量减少到 O(N / log₂N)，可以轻松节省对数级的执行时间并获得并行前缀问题的最优解。

以下是一个成本最优的部分和计算算法 parPrefix2：

const N = 8;
type Mass = array[1..N] of integer;
MassTemp = array[1..Floor(N / Log2(N))] of integer;
procedure parPrefix2(list: Mass; var prefixList: Mass);
    // 计算数组部分和的成本最优算法
    var i, j, k: integer;
    temp, prefixTemp: MassTemp;
begin
    for i := 1 to N do prefixList[i] := list[i];
    // 步骤 1
    parallel for i := 1 to Ceil(N / Log2(N)) do
        for j := 2 to Log2(N) do
        begin
            k := (i - 1) * Log2(N) + j;
            if (k <= N) then
                prefixList[k] := prefixList[k - 1] + list[k];
        end;
    // 步骤 2
    parallel for i := 1 to Floor(N / Log2(N)) do
        temp[i] := prefixList[i * Log2(N)];
    parPrefix(temp, prefixTemp);
    // 步骤 3
    parallel for i := 2 to Ceil(N / Log2(N)) do
        for j := (i - 1) * Log2(N) + 1 to i * Log2(N) do
            if (j <= N) then
                prefixList[j] := prefixTemp[i - 1] + prefixList[j];
end;

该算法各步骤的执行时间分析如下：
- 第一步和最后一步需要时间 O(log₂N)。
- 第二步计算数组 temp 的部分和，执行时间为 O(log₂(N / log₂N)) = O(log₂N)。

考虑到 p = O(N / log₂N)，可以得出并行前缀计算算法的第二个版本的成本估计为 Cp = O(N)。

并行算法：原理、模型与问题求解

4. 并行算法的性能分析与限制

在分析并行算法时，除了前面提到的加速比、效率和成本等概念，还需要考虑一些其他因素对算法性能的影响。

首先是 Amdahl 定律，它描述了由于存在无法并行化的顺序计算部分，对加速比的限制。设 f 是算法 A 中顺序计算的比例，那么在具有 p 个处理器的系统上使用算法 A 的加速比 Sp 满足不等式：
[Sp \leq \frac{1}{f + \frac{1 - f}{p}}]

证明过程如下：在多处理器上执行算法所需的时间由顺序操作 fT₁ 和可并行化的操作组成，即 (Tp = fT₁ + \frac{1 - f}{p}T₁)。因此，加速比的上限可以表示为：
[\left[Sp\right]_{max} = \frac{T₁}{fT₁ + \frac{(1 - f)T₁}{p}} = \frac{1}{f + \frac{1 - f}{p}}]

这个不等式表明，存在无法并行化的顺序计算会对加速比施加限制。即使使用超计算机（p → ∞），加速比也不会超过 (S_{\infty} = \frac{1}{f})。

不同 f 值下，加速比 Sp 与处理器数量 p 的关系如下图所示：
| f 值 | 加速比 Sp 与处理器数量 p 的关系 |
| — | — |
| f = 0.1 | 随着 p 增加，Sp 逐渐接近上限 10 |
| f = 0.05 | 随着 p 增加，Sp 逐渐接近上限 20 |
| f = 0.02 | 随着 p 增加，Sp 逐渐接近上限 50 |

从实际经验来看，对于广泛的计算任务，随着任务输入数据规模的增加，顺序计算的比例会减小。因此，在实践中，可以通过增加任务的计算复杂度来提高加速比。

然而，PRAM 模型虽然广泛应用，但也存在一些局限性。例如，它完全忽略了特定计算机系统架构中不同处理器的数据传输速率差异问题。

5. 算法结构的图形表示

算法解决问题的结构可以用有向图“操作 - 操作数”来图形化表示，这是一个无环有向图，记为 D = (V, E)，其中 V 是表示运行操作的顶点集合，E 是边的集合。边 (v_iv_j \in E) 当且仅当编号为 j 的操作使用编号为 i 的操作的结果。入度 (d^+(v_{k}^{(in)}) = 0)（k ∈ N）的顶点 (v_{k}^{(in)}) 用于数据输入操作，出度 (d^-(v_{l}^{(out)}) = 0)（l ∈ N）的顶点 (v_{l}^{(out)}) 对应输出操作。

graph LR
    classDef process fill:#E5F6FF,stroke:#73A6FF,stroke-width:2px
    A(v1<sup>(in)</sup>):::process --> B(v2):::process
    A --> C(v3):::process
    B --> D(v4):::process
    C --> D
    D --> E(v1<sup>(out)</sup>):::process

在假设使用 RAM 时，有向图 D 中每个顶点的计算时间为常数 τ 的情况下，超计算机上算法的执行时间 (T_{\infty}) 等于：
[T_{\infty} = \tau \times \max_{(i, j)} \left{\vert C(i, j) \vert \right}]
其中 (C(i, j)) 是从 (v_{i}^{(in)}) 到 (v_{j}^{(out)}) 的路径，最大值是对所有可能的 (i, j) 对取的。也就是说，时间 (T_{\infty}) 与连接某些 (v_{i}^{(in)}) 和 (v_{j}^{(out)}) 的最大路径中的顶点数量成正比。

6. 并行算法的实际应用与展望

并行算法在许多领域都有广泛的应用，如科学计算、数据处理、人工智能等。在科学计算中，并行算法可以加速复杂的数值模拟，如气象预报、分子动力学模拟等。在数据处理方面，并行算法可以提高大数据分析的效率，快速处理海量数据。在人工智能领域，并行算法可以加速深度学习模型的训练过程，提高模型的训练速度和性能。

然而，要实现高效的并行算法，还需要解决一些挑战。例如，如何更好地处理内存访问冲突，提高处理器之间的通信效率，以及如何设计更适合并行计算的算法结构等。未来，随着计算机硬件技术的不断发展，如多核处理器、GPU 等的广泛应用，并行算法的研究和应用将会更加深入和广泛。同时，新的并行计算模型和算法也将不断涌现，为解决各种复杂的计算问题提供更强大的工具。

总之，并行算法作为计算机科学领域的重要研究方向，具有巨大的发展潜力和应用前景。通过不断地研究和创新，我们可以充分发挥并行计算的优势，推动各个领域的发展和进步。