向量处理器（8）

4.3 稀疏矩阵

存在一些能够使得基于稀疏矩阵的程序在向量模式下运行的技术。在一个稀疏矩阵中，向量的元素通常是以紧凑的形式储存，并且以间接的方式被访问。我们会看到以下一个简化了的稀疏矩阵的代码：

                                                         do 100 i = 1, n

                     100                                    A(K(i)) = A(K(i)) + C(M(i))

这段代码利用 K 和 M 作为索引向量来给出 A 和 C 中的非零元素实现了数组 A 和 C 的稀疏向量求和和。（A 和 C 必须有相同数目个非零元素--本例中的 n）另外一个常见的稀疏矩阵的表示形式是用一个位向量来表示哪些位是非零元素以及一个稠密向量包含所有的非零元素。通常这两种表示形式会同时在一个程序里出现。在很多代码里都能看到稀疏矩阵的影子，并且根据不同的数据结构有很多中实现方法。

支持稀疏矩阵的最主要的机制是利用索引向量（index vector）的 scatter-gather 操作。这类操作的目标是支持在稠密表示（即没有非零元素）和正常表示（即包括非零元素）之间进行数据迁移。Gather 操作根据一个索引向量，通过把索引向量给出的偏移值加到基地址上来取出向量的元素。其输出是在一个在向量寄存器里的稀疏向量。在这些元素以稠密的形式被处理之后，稀疏向量要以扩展的形式通过 scatter 操作存回内存，使用同样的索引数组。对于这两个操作的硬件支持称之为 scatter-gather，并且在几乎所有的现在向量处理器上都能看到。VMIPS 提供了 LVI（load vector indexed）和 SVI（store vector indexed）指令来实现这两个操作。比如，假定 Ra，Rc，Rk 和 Rm 分别有前例中四个向量的起始地址，那么彼代码段的内层循环可以用以下的向量指令实现：

                                         LV                        Vk, Rk                       ;load K

                                         LVI                       Va, (Ra + Vk)          ;load A(K(I))

                                         LV                        Vm, Rm                    ;load M

                                         LVI                       Vc, (Rc + Vm)         ;load C(M(I))

                                         ADDV.D              Va, Va, Vc              ;add them

                                         SVI                       (Ra + Vk), Va          ;store A(K(I))

这个技术使得稀疏矩阵的代码可以以向量模式执行。简单的向量编译器不能自动向量化以上代码，因为编译器不会知道 K 中的元素的值互相是不相同的，因此不存在任何依赖关系 [1]。因此，程序员需要告诉编译器该循环可以以向量模式运行。

更为复杂的向量编译器可以自动向量化以上循环而不需要程序员的干预。这是通过插入运行时对于数据 dependency 的检查实现的。这种运行时检查是通过  Itanium 处理器中的 Advanced Load Address Table（ALAT）硬件机构的向量化软件版本实现的。ALAT 在 附录 G 中有描述。ALAT 硬件被一个软件哈希表所代替。该哈希表能够检测出在同一个 strip-mining 的迭代循环中的两个元素是否指向同一个地址。如果没有检测到 dependency，该 strip-mining 循环则可以以长度 MVL 来完成。如果检测到了 dependency，向量的长度则被重置为较小的一个可以避免 dependency 的值，而留下剩下的部分给下一个循环执行。虽然这种机制给执行循环增加了很多软件开销，但是这种开销还是会被更为常见的没有依赖的情况所均摊，因此这个循环仍然会比标量代码快得多（当然会比程序员直接指出可以向量化的情况慢得多）。

近来大多数的超级计算机都有 scatter-gather 的能力。这种操作比有跨度的操作更为之慢因为实现起来更复杂，而且更容易出现 bank 冲突，但是比之标量版本，则要快很多。如果一个矩阵的稀疏程度改变了，必须重新计算索引向量。很多处理器提供了快速计算所以向量的方法。VMIPS 中的 CVI 指令可以根据一个给定的跨度值（m）来创建一个所以向量。其各个元素值为0，m，2×m，...，63×m。一些处理器提供一条创建压缩形式的索引向量的指令。该向量中各元素值对应于掩码寄存器中相应位置为 1 的元素。另外一些则提供压缩向量的指令。在 VMIPS 中，我们定义 CVI 指令为总是根据向量掩码来创建一个压缩过的索引向量。当所有的掩码都为 1 的时候，则创建一个标准的索引向量。

索引化的 L/S 操作以及 CVI 指令提供了支持条件向量执行的一种新方法。以下是我们在上一小节中循环的另一种实现：

                                         LV                        V1, Ra                     ;load vector A into V1

                                         L.D                       F0, #0                      ;load FP zero into F0

                                         SNEVS.D             V1, F0                    ;sets the VM to 1 if V1(i) != F0

                                         CVI                      V2, #8                      ;generates indices in V2

                                         POP                     R1, VM                   ;find the number of 1's in VM

                                         MTC1                  VLR, R1                   ;load vector-length register

                                         CVM                                                      ;clears the mask

                                         LVI                       V3, (Ra + V2)         ;load the nonzero A elements

                                         LVI                       V4, (Rb + V2)         ;load corresponding B elements

                                         SUBV.D               V3, V3, V4             ;do the subtract

                                         SVI                       (Ra + V2), V3         ;store A back

至于究竟是使用 scatter-gather 的版本更好还是使用条件执行的版本更好取决于该条件测试满足的频率以及这些操作的开销。不考虑 chaining 的话，第一个版本（前一小节）的耗时是 5n + c1。第二个版本，也即采用每一个周期能执行对于一个元素索引化 L/S 的版本，的执行时间是4n + 4fn + c2 [2]，其中 f 是条件测试满足（也即A(i) != 0）的比率。如果我们假设 c1 和 c2 差不多，或者说他们都远小于 n，我们可以求出什么时候第二个版本会更好。

                                         Time1 = 5n

                                         Time2 = 4n + 4fn

如果我们要 Time1 > Time2，那么

                                         5n > 4n + 4fn

                                         1/4 > f

也就是说，如果小于 1/4 的元素是非零元素，那么第二个版本更好。在很多情况下，条件满足的比率要小得多。如果可以索引向量可以被重用，或者在 if 语句下的向量语句 [3] 的数目增加，scatter-gather 的优势会显著增加。

---------------大家好，我是分割线---------------

[1] 如果 K 中有两个元素的值相等，那么它们会指向对于 A 中同一个元素进行操作。简单的编译器会认为这两个操作之间有相互数据依赖而不能自动并行化。

[2] 为什么是 4n + 4fn？在第二个版本中，需要完全执行（也即需要遍历向量中所有的元素）的向量指令是 LV，SNEVS.D，CVI，POP；需要部分执行（也即只需要对满足条件的向量元素进行操作）的指令是 2 条 LVI，SUBV.D，SVI。条件满足的比率是 f。

[3] 也就是满足条件测试后所需要执行的那部分语句。拿本例来说，如果 A(i) != 0 的情况下不仅仅是做简单的减法而是执行一系列复杂的操作，那么索引化L/S 的优势更明显。原因在于索引化 L/S 相比原来只用掩码的版本的优点正在于减少了 if 语句下的那些指令的执行时间（即只对满足条件的元素执行操作）。如果 if 语句下的那些指令数目增大，这个优势会被扩大。读者可以用上面给出的分析模型定性分析。