本章节介绍 the minimum degree algorithm ，我个人将其翻译为最小度数算法，我们这里看消元图的理论基础。

在这一章节中我们研究的排序算法称为 minimum degree algorithm (Rose) ，这是一个启发式算法，用于对 \(\mathbf{A}\) 排序减少其分解时产生的 fill-in 数目。

Cholesky 分解

设 \(\mathbf{A}\) 是一个对称稀疏矩阵，其中的非零结构 (nonzero structure) 定义为 \[ Nonz(\mathbf{A}) = \{ \{ i,j \} \mid a_{ij} \neq 0 \text{ and } i \neq j \} \] 假设 \(\mathbf{A}\) 分解为 \(\mathbf{LL^{T}}\) ，而 filled matrix \(\mathbf{F(A)＝L+L^{T}}\)，下文中我们用 \(\mathbf{F}\) 替换 \(\mathbf{F(A)}\) ，其相应的非零结构为 \[ Nonz(\mathbf{F}) = \{ \{ i,j \} \mid l_{ij} \neq 0 \text{ and } i \neq j \} \] 在本文中，我们假设精确的数值消除并不存在（就是说任意两个浮点数相减不能得到0），因此对于一个给定的非零结构 \(Nonz(\mathbf{A})\)，其相应的 \(Nonz(\mathbf{F})\) 的结构就能完全确定。

这个假设立刻推导出 \[ Nonz(\mathbf{A}) \subset Nonz(\mathbf{F}) \] 矩阵 \(\mathbf{A}\) 的 fill 可以定义为 \[ Fill(\mathbf{A}) = Nonz(\mathbf{F}) - Nonz(\mathbf{A}) \] 举个例子，假设存在下图中的矩阵，其中的 fill-in 用加号表示，其相应的集合为 $$ \begin{aligned} \operatorname{Nonz}(\boldsymbol{A}) &=\{\{1,5\},\{1,8\},\{2,4\},\{2,5\},\{3,8\},\{4,7\},\{5,6\},\{6,8\},\{8,9\}\} \\ \operatorname{Fill}(\boldsymbol{A}) &=\{\{4,5\},\{5,7\},\{5,8\},\{6,7\},\{7,8\}\} . \end{aligned} $$

Elimation Graph Model

我们现在联系 \(\mathbf{A}\) 的高斯消元，与相应的图 \(\mathcal{G}^{\mathbf{A}}\) 的变化。我们回顾一下分解的 outer product form ，第一步为 $$ \begin{aligned} \boldsymbol{A} &=\boldsymbol{A}_{0}=\boldsymbol{H}_{0}=\left(\begin{array}{cc} d_{1} & \boldsymbol{v}_{1}^{T} \\ \boldsymbol{v}_{1} & \overline{\boldsymbol{H}}_{1} \end{array}\right) \\ &=\left(\begin{array}{cc} \sqrt{d_{1}} & \mathbf{0} \\ \frac{\boldsymbol{v}_{1}}{\sqrt{d_{1}}} & \boldsymbol{I}_{n-1} \end{array}\right)\left(\begin{array}{cc} 1 & \mathbf{0} \\ \mathbf{0} & \boldsymbol{H}_{1} \end{array}\right)\left(\begin{array}{cc} \sqrt{d_{1}} & \frac{\boldsymbol{v}_{1}^{T}}{\sqrt{d_{1}}} \\ \mathbf{0} & \boldsymbol{I}_{n-1} \end{array}\right) \\ &=\boldsymbol{L}_{1} \boldsymbol{A}_{1} \boldsymbol{L}_{1}^{T}, \end{aligned} $$ 其中 $$ \boldsymbol{H}_{1}=\overline{\boldsymbol{H}}_{1}-\frac{\boldsymbol{v}_{1} \boldsymbol{v}_{1}^{T}}{d_{1}} $$ 这一步会递归地应用于 $\mathbf{H}_{1}$， $\mathbf{H}_{2}$ 等。由于我们假设不存在精准的消除，因此根据上面的公式，如果 \(\mathbf{\bar{H}_{1}}\) 中的 \(\{j,k\}\) 元素不为0，或者 \((v_{1})_{j} \neq 0 \text{ and } (v_{1})_{k} \neq 0\) ，\(\mathbf{H}_{1}\) 矩阵的 \(\{j,k\}\) 元素才不为0。当第一个条件不成立，而第二个条件成立时，便产生了 fill-in 。下图中展示了这个情况

　

当第一步分解完成之后，我们对 \(\mathbf{H}_{1}\) 矩阵继续进行分解。

现在我们看从 \(\mathbf{H}_{0}\) 到 \(\mathbf{H}_{1}\) 的转变过程中，其相应的图从 \(\mathcal{G}^{\mathbf{H}_{0}}\) 到 \(\mathcal{G}^{\mathbf{H}_{1}}\)的变化，下图是一个例子。生成图 \(\mathcal{G}^{\mathbf{H}_{1}}\)的过程如下

1. 删除节点 \(x_{1}\) 及其相应的边
2. 添加边，使得图 \(\mathcal{G}^{\mathbf{H}_{1}}\)中的 \(Adj(x_{1})\) 彼此之间相邻

　

因此，就像 Rose 观察到的，对称矩阵的高斯消元（上面的 \(\mathbf{H}_{1}\)， \(\mathbf{H}_{2}\) 等其实也一样就是高斯消元的结果 ）可以理解为生成一系列 elimination graphs \[ \mathcal{G}_{i}^{\alpha}=\mathcal{G}^{\boldsymbol{H}_{i}}=\left(X_{i}^{\alpha}, E_{i}^{\alpha}\right), \quad i=1,2, \cdots, n-1 \] 其中 \(\mathcal{G}_{i}^{\alpha}\) 是从 \(\mathcal{G}_{i-1}^{\alpha}\) 中得到的，过程就和上面描述的一样。这里 \(\alpha\) 是所有节点的标签，当我们已知所有节点的标签时，我们用标记 \(\mathcal{G}_{i}\) 替换 \(\mathcal{G}_{i}^{\alpha}\) 。在上图的消元过程，比如从 \(\mathcal{G}_{1}\) 消除 \(x_{2}\) 的过程中，产生了三个 fill-in 的边 \(\{x_{3}, x_{4}\}, \{x_{3}, x_{6}\}, \{x_{4}, x_{6}\}\) 。

设 \(\mathbf{L}\) 是 \(\mathbf{A}\) 的分解因子，定义 \(\mathcal{G}^{\mathbf{A}}\)的 filled graph 为 \(\mathcal{G}^{\mathbf{F}} = (X^{\mathbf{F}},E^{\mathbf{F}})\) ，其中 \(\mathbf{F = L + L^{T}}\) 。这里所有边的集合 \(E^{\mathbf{F}}\)包含了 \(E^{\mathbf{A}}\) 中的所有边，加上分解过程中新增的边。

\(E^{\mathbf{F}}\) 和 \(E^{\mathbf{A}}\) 的关系存在以下引理 (Parter)

引理：当且仅当 \(\{ x_{i}, x_{j}\} \in E^{\mathbf{A}}\)，或者 \(\{ x_{i}, x_{k}\} \in E^{\mathbf{F}}\) and \(\{ x_{k}, x_{j}\} \in E^{\mathbf{F}}\) for some \(k < \min\{i,j\}\) ，才满足 \(\{ x_{i}, x_{j}\} \in E^{\mathbf{F}}\)

根据这个引理，上图中得到的 \(\mathcal{G}^{\mathbf{F}}\) 如下。发现了 \(\mathcal{G}^{\mathbf{F}}\) ，我们就可以得到 \(\mathbf{L}\) 矩阵的结构。

　

Modelling Elimination By Reachable Sets

上面的章节中定义了 elimination graphs 的序列 \[ \mathcal{G}_{0} \rightarrow \mathcal{G}_{1} \rightarrow \cdots \rightarrow \mathcal{G}_{n-1} \] 然后提供了一种递归的方式来确认 \(E^{\mathbf{F}}\) ，我们这一章节的目标是采用 reachable sets 的概念来确认这些步骤。

让我们首先研究在上图中 fill edge \(\{ x_{4}, x_{6} \}\) 的形成过程。在 \(\mathcal{G}_{1}\) ，存在路径： \[ (x_{4}, x_{2}, x_{6}) \] 因此当 \(x_{2}\) 被消除，被创建了 \(\{ x_{4}, x_{6} \}\) 这条边。然而，\(\{ x_{2}, x_{6} \}\) 并不在原始的图中，这条边是由于路径 \((x_{2}, x_{1}, x_{6})\) 中的 \(x_{1}\) 从 \(\mathcal{G}_{0}\) 中消除而创建得到的。

联合这两个，我们得到 \(\mathcal{G}_{0}\) 中的路径 \((x_{4}, x_{2}, x_{1}, x_{6})\) 才是生成 fill edge \(\{ x_{4}, x_{6} \}\) 的真正原因，因此我们引入 reachable sets 的概念。

设 \(S\) 是节点集合的一个真子集，\(x \notin S\) 。如果存在一个从 y 到 x 的路径 \((y, v_{1}, \cdots, v_{k},x)\) ，其中 \(v_{i} \in S, 1 \leq i \leq k, k \geq 0\)，则我们称节点 \(x\) 是 be reachable from a node \(y\) through \(S\) ，注意由于 \(k\) 可以为 0，因此 \(y\) 的任何不在 \(S\) 中的相邻节点均 be reachable from \(y\) through \(S\)

这些节点的集合就是 reachable set of y through S ，标记为 \(Reach(y, S)\) ，定义为 \[ \operatorname{Reach}(y, S)=\{x \notin S \mid x \text { is reachable from } y \text { through } S\} . \] 为了说明 reachable set 的概念，以下图为例，如果 \(S = \{s_{1}, s_{2},s_{3},s_{4}\}\) ，我们有 \[ \operatorname{Reach}(y, S)=\{ a,b,c\} \] 　因为存在以下路径 \(\{y,s_{2},s_{4},a\}\)，\(\{y,b\}\)，\(\{y,s_{1},c\}\)

下面的定理说明通过 reachable set 来创建的 filled graph

定理： \[ E^{\boldsymbol{F}}=\left\{\left\{x_{i}, x_{j}\right\} \mid x_{j} \in \operatorname{Reach}\left(x_{i},\left\{x_{1}, x_{2}, \cdots, x_{i-1}\right\}\right)\right\} . \] 证明：假设 \(x_{j} \in \operatorname{Reach}\left(x_{i},\left\{x_{1}, x_{2}, \cdots, x_{i-1}\right\}\right)\) ，根据定义，在 \(\mathcal{G}_{0}\) 中存在一条路径 \((x_{i},y_{1},\cdots, y_{t},x_{j})\) ，其中 \(y_{k} \in \left\{x_{1}, x_{2}, \cdots, x_{i-1}\right\}\) 。如果 \(t = 0\) 或者 \(t = 1\) ，根据上面的引理就可以直接得到 \(\left\{x_{i}, x_{j}\right\} \in E^{\mathbf{F}}\)。如果 \(t > 1\) ，我们可以从 \(y_{1},\cdots, y_{t}\) 中逐渐消除节点，使得 \(t = 0\) ，从而得到 \(\left\{x_{i}, x_{j}\right\} \in E^{\mathbf{F}}\)。

反过来，假设 \(\left\{x_{i}, x_{j}\right\} \in E^{\mathbf{F}}\) ，并且 \(i \leq j\) 。根据前面的推理，则 \(\{ x_{i}, x_{j}\} \in E^{\mathbf{A}}\)，或者 \(\{ x_{i}, x_{k}\} \in E^{\mathbf{F}}\) and \(\{ x_{k}, x_{j}\} \in E^{\mathbf{F}}\) for some \(k < \min\{i,j\}\) 。如果是满足第一个条件 \(\{ x_{i}, x_{j}\} \in E^{\mathbf{A}}\)，则 \(x_{j} \in \operatorname{Reach}\left(x_{i},\left\{x_{1}, x_{2}, \cdots, x_{i-1}\right\}\right)\) 成立；如果是满足第二个条件，即 \(\{ x_{i}, x_{k}\} \in E^{\mathbf{F}}\) and \(\{ x_{k}, x_{j}\} \in E^{\mathbf{F}}\) for some \(k < \min\{i,j\}\) ，这说明在 \(\mathcal{G}_{0}\) 中 \(x_{i}\) 和 \(x_{k}\)， \(x_{k}\) 和 \(x_{j}\) 均存在一条路径连接，二者可以连通起来，说明 存在某条通过 \(x_{k}\) 连接 \(x_{i}\) 和 \(x_{j}\) 的路径 \((x_{i},x_{k},x_{j})\) ，因此 \(x_{j} \in \operatorname{Reach}\left(x_{i},\left\{x_{1}, x_{2}, \cdots, x_{i-1}\right\}\right)\) 成立。

证明完毕。

采用矩阵的术语，\(\operatorname{Reach}\left(x_{i},\left\{x_{1}, x_{2}, \cdots, x_{i-1}\right\}\right)\) 就是在列向量 \(\mathbf{L}_{*i}\) 中的非零元素的行号。举个例子，假设上图 5.2.5 中的图排序如下

如果 \(S_{i} = \left\{x_{1}, x_{2}, \cdots, x_{i-1}\right\}\) ，根据上面的定理，我们可以轻松得到下面的 reachable set \[ \begin{aligned} \operatorname{Reach}\left(x_{1}, S_{0}\right) &=\left\{x_{5}, x_{8}\right\} \\ \operatorname{Reach}\left(x_{2}, S_{1}\right) &=\left\{x_{4}, x_{5}\right\} \\ \operatorname{Reach}\left(x_{3}, S_{2}\right) &=\left\{x_{8}\right\} \\ \operatorname{Reach}\left(x_{4}, S_{3}\right) &=\left\{x_{5}, x_{7}\right\} \\ \operatorname{Reach}\left(x_{5}, S_{4}\right) &=\left\{x_{6}, x_{7}, x_{8}\right\} \\ \operatorname{Reach}\left(x_{6}, S_{5}\right) &=\left\{x_{7}, x_{8}\right\} \\ \operatorname{Reach}\left(x_{7}, S_{6}\right) &=\left\{x_{8}\right\} \\ \operatorname{Reach}\left(x_{8}, S_{7}\right) &=\left\{x_{9}\right\} \\ \operatorname{Reach}\left(x_{9}, S_{8}\right) &=\phi . \end{aligned} \] 因此我们可以得到下面的 \(\mathbf{L}\) 矩阵

因此我们可以通过 \(\mathbf{A}\) 的结构而直接得到 \(\mathbf{L}\) 矩阵的结构（这种方式比 envelope region 的方式更加精确，fill-in 更少）。更重要的是，我们可以更方便的得到 elimination graphs 。设 \(\mathcal{G}_{0} \rightarrow \mathcal{G}_{1} \rightarrow \cdots \rightarrow \mathcal{G}_{n-1}\) 为一组 elimination graphs ，我们有：

定理：设 \(y\) 是一个 elimination graph \(\mathcal{G}_{i} = (X_{i},E_{i})\) 中的一个节点，邻近于 \(\mathcal{G}_{i}\) 中的节点 \(y\) 的节点集合为： \[ \operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{i}\right\}\right) \] 其中，这里的 \(\operatorname{Reach}\) 操作符是用于原始的图 \(\mathcal{G}_{0}\) 。

证明：证明可以通过对 \(i\) 采用归纳法证明得到，设 \(i = 1\) ，易证明成立；假设 \(i =k\) 时成立，即邻近于 \(\mathcal{G}_{k}\) 中的节点 \(y\) 的节点集合为 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right)\) ；当 \(i = k+1\) 时， \(\mathcal{G}_{k+1}\) 进一步剔除了 \(x_{k+1}\) ，我们需要证明邻近于 \(\mathcal{G}_{k+1}\) 中的节点 \(y\) 的节点集合为 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k+1}\right\}\right)\) 。此时如果 \(y \notin \operatorname{Reach}\left(x_{k+1},\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right)\) , 也就是说在 \(\mathcal{G}_{k}\) 中 \(y\) 和 \(x_{k+1}\) 不相邻，此时 \(\mathcal{G}_{k+1}\) 中的节点 \(y\) 的相邻列表不变，仍为 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right)\) ，也等于 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k+1}\right\}\right)\) 。如果 \(y \in \operatorname{Reach}\left(x_{k+1},\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right)\) ，此时 \(y\) 的相邻节点会减去 \(x_{k+1}\) , 新增 \(x_{k+1}\) 在 \(\mathcal{G}_{k}\) 中 除 \(y\) 之外的相邻节点，即此时 \(\mathcal{G}_{k+1}\) 中的节点 \(y\) 的相邻列表为 \((\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right) - \{x_{k+1}\}) \cup (\operatorname{Reach}\left(x_{k+1},\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right) - \{y\})\) ，易得其中第一部分为 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k}, x_{k+1}\right\}\right)\) 中路径中不包含 \(x_{k+1}\) 的部分 ，对于第二部分，设 \(t \in (\operatorname{Reach}\left(x_{k+1},\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right) - \{y\})\) , 由于 \(y \in \operatorname{Reach}\left(x_{k+1},\left\{x_{1}, x_{2}, \cdots, x_{k}\right\}\right)\) ，易得 \(t \in \operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k}, x_{k+1}\right\}\right)\) ，因此第二部分为 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k}, x_{k+1}\right\}\right)\) 路径中包含 \(x_{k+1}\) 的部分 ，因此第一部分和第二部分的并集就正好为 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k}, x_{k+1}\right\}\right)\) ，因此得证 邻近于 \(\mathcal{G}_{k+1}\) 中的节点 \(y\) 的节点集合为 \(\operatorname{Reach}\left(y,\left\{x_{1}, x_{2}, \cdots, x_{k+1}\right\}\right)\) 。

证明完毕。

我们用图 5.2.3 来重新举一个例子，我们考虑 图 \(\mathcal{G}_{0}\) 和 图 \(\mathcal{G}_{2}\) 如下：

设 \(S_{2} = {x_{1}, x_{2}}\) ，我们易得： $$ \[\begin{aligned} \operatorname{Reach}\left(x_{3}, S_{2}\right) &=\left\{x_{4}, x_{5}, x_{6}\right\} \\ \operatorname{Reach}\left(x_{4}, S_{2}\right) &=\left\{x_{3}, x_{6}\right\} \\ \operatorname{Reach}\left(x_{5}, S_{2}\right) &=\left\{x_{3}, x_{6}\right\} \\ \operatorname{Reach}\left(x_{6}, S_{2}\right) &=\left\{x_{3}, x_{4}, x_{5}\right\} \\ \end{aligned}\]

$$ 因此我们就得到 图 \(\mathcal{G}_{2}\) 的结构。

我们采用 reachable sets 得到的 图 \(\mathcal{G}_{i}\) 称为 implicit model for elimation；而上一节中消除节点 \(x_{1}\) 及其相应的边 添加边，使得图 \(\mathcal{G}^{\mathbf{H}_{1}}\)中的 \(Adj(x_{1})\) 彼此之间相邻，这种方式称为 explicit model for elimation。

参考文献

1. George A, Liu J, Ng E. Computer solution of sparse linear systems[J]. Oak Ridge National Laboratory, 1994.