稀疏矩阵

Sparse Matrix
大部分元素都是零的矩阵

在科学计算、数据分析、机器学习等领域，稀疏矩阵经常出现，因为很多实际问题中的矩阵表示会有大量的零元素。

处理稀疏矩阵时，通常不会直接存储所有的零值，而是只存储和操作非零元素，这样可以大大减少内存使用和计算时间。

稀疏矩阵在处理大规模数据集时非常有用，因为它们可以显著减少内存占用和计算资源，从而提高算法的效率。

稀疏矩阵几乎产生于所有的大型科学工程计算领域，包括计算流体力学、统计物理、电路模拟、图像处理、纳米材料计算等。

存储方式

压缩稀疏行（Compressed Sparse Row, CSR）：将矩阵的行表示为三个数组：一个存储非零值，一个存储每行的非零值的列索引，一个存储每行的起始位置。这种方法适合于行操作和矩阵运算。

压缩稀疏列（Compressed Sparse Column, CSC）：与CSR类似，但是是按列来存储数据，适合于列操作。

块压缩稀疏行（Block Compressed Row Storage, BCRS）：是CSR的一种变体，将矩阵分割成块，每个块作为一个单元存储。

带宽与带状矩阵

许多稀疏矩阵的非零元素集中在主对角线附近。若存在半带宽 $w$ 使得

a_{i j} = 0 (| i - j | > w),

则称矩阵具有带状结构。对角矩阵可看作最窄的带状矩阵，三对角矩阵只在主对角线及其上下相邻对角线上有非零元素。

带状结构重要，是因为高斯消元可以只处理带内元素。若消元不破坏带状结构，前向消元的主要工作量可由稠密矩阵的

\frac{1}{3} n^{3}

降到约

w^{2} n,

而一次前代加回代约为 $2 w n$ 。当 $w$ 不随 $n$ 增大时，大规模带状系统可以非常便宜。

填充

稀疏矩阵的困难在于，原来的零元素不一定在消元后继续为零。消元过程中由零变成非零的位置称为填充。填充会增加 $L, U$ 或 Cholesky 因子的非零数，使存储和运算都上升。

一个 $3 \times 3$ 例子可以说明填充：

A = [\begin{matrix} 1 & 1 & 1 \\ - 2 & 1 & 0 \\ - 2 & 0 & 2 \end{matrix}] .

用第一行消去下方第一列后，第三行第二列会从 $0$ 变成 $2$ ，即

a_{32} = 0 ⟶ a_{32} = 2.

这个新非零不是输入矩阵原有的结构，而是消元顺序造成的。稀疏直接法的核心目标之一就是减少这种新非零。

图解释与重排序

稀疏矩阵常用图来理解：把每个未知量或方程看作节点，若 $a_{i j} \neq 0$ ，就在节点 $i, j$ 之间建立连接。消元一个节点时，它的邻居之间可能需要补上新的边；这些新边对应因子中的填充。

因此，稀疏求解器常在数值分解前先做重排序。设 $P$ 为置换矩阵，重排序后的系统可写成

P^{T} A P 或 A P .

目标不是改变问题的解，而是改变消元顺序，使填充更少、带宽更窄或前沿矩阵更小。常见思想包括：

最小度：优先消去当前非零邻居最少的节点。
AMD：用近似最小度降低精确最小度的搜索成本。
colamd：面向非对称问题的列近似最小度排序，选择列置换 $P$ 来减少 $A P$ 在稀疏 LU 中的填充。
symamd：面向对称问题的近似最小度排序，选择同一个置换同时作用在行和列上，以减少 $P^{T} A P$ 在 Cholesky 或对称分解中的填充。
带宽缩减：让相关节点在编号上更接近，适合明显局部连接的问题。

两者的差别在于结构假设不同：对称重排序必须同时重排行和列，才能保持矩阵图的无向结构与对称性；非对称列重排序只改变列顺序，主要服务于一般稀疏 LU 的列主元结构和填充控制。

寻找绝对最小填充通常代价过高，所以实际软件使用启发式排序，在质量和预处理成本之间折中。

稀疏分解

稠密矩阵的 LU 或 Cholesky 分解主要关心浮点运算；稀疏分解还必须先关心非零模式。典型流程包括：

符号分析：根据非零位置和排序预估因子的结构、内存和填充。
数值分解：在已知结构上计算实际数值，必要时结合主元策略。
三角求解：利用稀疏 $L, U$ 或 Cholesky 因子求解一个或多个右端项。

对称正定矩阵常用稀疏 Cholesky：

P A P^{T} = L L^{T} .

一般非对称矩阵常用稀疏 LU：

P A Q = L U .

其中 $P, Q$ 记录行列置换。稀疏 LU 还要在数值稳定性和填充之间折中：选择较大主元有利于控制舍入误差，但某些行交换可能制造更多填充。

与逆矩阵的区别

稀疏矩阵本身稀疏，并不意味着 $A^{- 1}$ 也稀疏。三对角矩阵的逆往往是满矩阵。若为了求解 $A x = b$ 而显式形成 $A^{- 1}$ ，可能把一个只需 $O (n)$ 存储和计算的问题变成 $O (n^{2})$ 存储和更高成本的问题。

例如下面的 $4 \times 4$ 三对角矩阵只有 10 个非零元，但它的逆矩阵没有零元：

A = [\begin{matrix} 1 & - 1 & 0 & 0 \\ - 1 & 2 & - 1 & 0 \\ 0 & - 1 & 2 & - 1 \\ 0 & 0 & - 1 & 2 \end{matrix}], A^{- 1} = [\begin{matrix} 4 & 3 & 2 & 1 \\ 3 & 3 & 2 & 1 \\ 2 & 2 & 2 & 1 \\ 1 & 1 & 1 & 1 \end{matrix}] .

这说明稀疏性通常应保存在矩阵和分解因子中，而不是寄希望于逆矩阵继续稀疏。

因此，在稀疏问题中通常保存因子并做三角求解：

L y = b, U x = y,

而不是计算 $x = A^{- 1} b$ 。这既保留结构，也减少舍入和数据移动带来的额外风险。

例子

一维二阶差分离散常产生三对角矩阵：

[\begin{matrix} 2 & - 1 & 0 & 0 \\ - 1 & 2 & - 1 & 0 \\ 0 & - 1 & 2 & - 1 \\ 0 & 0 & - 1 & 2 \end{matrix}] .

按自然顺序消元时，它的因子仍接近双对角或窄带结构。若随机打乱行列顺序，同一个问题的矩阵仍然稀疏，但带宽可能变大，消元中会产生更多填充。这说明稀疏求解的速度不只取决于非零个数，还取决于非零位置。

二维网格问题更典型。每个网格点只连接上下左右邻居，原矩阵很稀疏；但无论按行编号还是按列编号，都可能在消元后产生远距离填充，所以需要更细致的排序策略。

边界条件

稀疏矩阵的优势有明确边界。若非零元素分布接近随机，填充可能很快让因子接近稠密；若矩阵很小，稀疏格式的索引开销可能超过收益；若操作需要大量随机访问，CSR 或 CSC 的效率会依赖访问方向；若数值主元很小，稳定性要求可能迫使额外行交换，从而牺牲一部分稀疏结构。

高斯消元 - 稀疏矩阵求解中最需要控制填充和主元的直接法基础。
部分主元法 - 稀疏 LU 中用于稳定性的行交换策略，但可能增加填充。
舍入误差 - 小主元和大乘子会让稀疏分解中的数值误差放大。
Cholesky分解 - 对称正定稀疏系统常用的分解形式。
图论 - 用节点和边描述非零结构、填充和重排序。