理解稀疏 Cholesky 消去树

在数值线性代数领域，Cholesky 分解是求解对称正定系统的基石。然而，在处理稀疏矩阵（即大多数元素为零）时，应用稠密分解算法在计算上是浪费的，且极其消耗内存。挑战在于“填充”（fill-in）：即原始矩阵 $A$ 中的零元素在生成的下三角矩阵 $L$ 中变为非零元素的现象。

为了管理这一点，工程师和数学家使用消去树（Elimination Tree）。这种结构允许我们精确预测填充会发生在哪里，并确定最佳的操作顺序，将复杂的任务依赖图转换为易于管理的树结构。本文将探讨如何直接从右向 Cholesky 算法中推导出消去树。

Cholesky 的机制与填充问题

要理解消去树，我们必须首先观察稠密右向 Cholesky 算法。对于每个主元 $k$，该过程涉及三个主要步骤：

  1. 主元分解：计算对角线元素 $L[k][k] = \sqrt{A[k][k]}$。
  2. 列缩放：缩放主元下方的列。
  3. 秩-1 更新：更新剩余矩阵：$L[i][j] -= L[i][k] * L[j][k]$。

正是这第三步——秩-1 更新——引入了填充。具体来说，如果我们有非零条目 $L[i][k]$ 和 $L[j][k]$（其中 $k < j \le i$），该更新必然会在 $L[i][j]$ 处创建一个或修改一个非零条目。

从任务 DAG 到消去树

如果我们映射出稀疏矩阵所需的每一个操作，最初会得到一个任务依赖的有向无环图（DAG）。虽然这个 DAG 告诉了我们所需的一切信息，但它通常是冗余的。

由于上述结构规则（$L[i][k] \ \neq 0$ 且 $L[j][k] \ \neq 0 \implies L[i][j] \ \neq 0$），依赖图中的许多边是隐含的。例如，如果第 0 列依赖于第 1 列，且第 1 列依赖于第 2 列，那么第 0 列对第 2 列的依赖关系已经包含在内了。

通过移除这些冗余边，DAG 会坍缩成一棵列消去树（Column Elimination Tree）。这棵树是一个 $O(n)$ 的数据结构，提供了两个关键信息：

  1. 填充预测：精确预测 $L$ 中哪些条目将是非零的。
  2. 任务调度：必须处理列的精确顺序。

实现消去树

符号分解

在进行数值计算之前，我们会进行“符号分解”。这一步使用消去树为 $L$ 预分配内存，而无需实际计算数值。算法遍历 $A$ 的行，并沿着消去树向上行走以标记所有祖先，从而确保在 L_col 结构中考虑到每一个潜在的填充位置。

数值分解

有了预先计算好的非零模式，数值分解就变成了遍历预分配索引的问题。秩-1 更新仅在已知的非零条目上执行，从而避免了在计算的内层循环中检查零值的开销。

计算树结构

定义消去树非常简单：列 $k$ 的父节点是满足 $L[j,k] \neq 0$ 的最小索引 $j > k$。

为了仅使用 $A$ 的初始非零元素来高效地计算此树，我们可以使用增量法。通过按递增顺序处理行 $r$，并维护一个 ancestor 数组，我们可以追踪从列 $c$ 到当前行 $r$ 的路径。如果路径终止，当前行 $r$ 成为该路径中最后一个节点的父节点。这确保了父节点是第一个需要该列作为祖先的后续行。

结论

通过将消去树建立在右向 Cholesky 算法而非抽象图论的基础上，线性代数与软件实现之间的关系变得清晰可见。消去树不仅是一个理论构造，更是一个实用的工具，它将稠密分解的二次复杂度转化为针对输入矩阵稀疏性定制的高效过程。

Sources