1. 基本概念

1.1 矩阵的图表示

对于  对称矩阵 ，定义其无向图 ：

• 顶点集  对应矩阵的行/列索引

• 边集 

由于 对称， 。

例子：

对应的图：边集 ，形成一条路径 1-2-3-4 加上 1-4 边，构成一个4-cycle加一条对角线。

2. Cholesky分解与填充图

2.1 填充现象

当对稀疏矩阵  进行 Cholesky 分解  时， 中可能出现一些  中原本为零的位置变成非零。这些新出现的非零元称为填充（fill-in）。

填充图   定义为Cholesky因子 的对称模式对应的无向图：

• 顶点集不变

• 边集包含所有原始边以及所有填充边

2.2 填充产生的图论解释

定理：在 Cholesky 分解过程中，当消去顶点  时，会在其所有邻居之间产生新的边（填充）。即如果  和  都是  的邻居，且 ，则边  会出现在填充图中。

这个过程的图论解释是：消元过程将图逐步转化为弦图。

3. 弦图与填充图的关系

3.1 弦图定义

        一个无向图称为弦图（chordal graph），如果它不包含长度≥4的无弦环。换句话说，任何长度  的环都至少有一条对角线（弦）。

3.2 重要定理

定理：Cholesky 分解的填充图  是一个弦图。而且，它是包含原始图  的最小弦图（即弦图的完备化）。

定理：原始图  是弦图当且仅当存在一个排列（消元顺序）使得 Cholesky 分解不产生任何填充。

4. 消去树

4.1 定义

消去树（elimination tree）是描述 Cholesky 分解中依赖关系的有向树（或森林）：

• 顶点集与矩阵相同

• 对每个 ，如果  且对所有 ，，则  是   的父节点

• 更简单的定义：顶点jj的父节点是 

4.2 性质

1. 消去树反映了矩阵分解的数据依赖关系

2. 如果 ，则  是  的祖先（不一定是父节点）

3. 消去树可以用于指导并行计算

5. 在Cholesky分解中的应用

5.1 符号分解

在实际计算前，先进行符号分解，预测填充位置：

def symbolic_factorization(G, n):
    """
    符号分解：预测Cholesky因子的非零模式
    G: 原始图的邻接表
    n: 顶点数
    返回：L的非零模式（邻接表）
    """
    # 初始化L的非零模式为原始图的下三角部分
    L_pattern = [set() for _ in range(n)]
    for i in range(n):
        for j in G[i]:
            if j < i:
                L_pattern[i].add(j)
    
    # 模拟消元过程
    reach = [set() for _ in range(n)]  # 可达集
    for k in range(n):
        # 计算顶点k的可达集
        reach_set = set()
        for j in L_pattern[k]:
            reach_set.update(reach[j])
        reach_set.update(L_pattern[k])
        reach[k] = reach_set
        
        # 填充产生：k的邻居之间产生新边
        neighbors = [j for j in L_pattern[k]] + [k]
        for i in neighbors:
            for j in neighbors:
                if i > j and j not in L_pattern[i]:
                    L_pattern[i].add(j)
    
    return L_pattern

 

5.2 重排序

通过重新排列矩阵的行列（图的重标号），可以减少填充：

最小度算法（Minimum Degree Algorithm）：

def minimum_degree_ordering(G, n):
    """
    最小度排序：每次选择度最小的顶点消去
    返回：新的顶点顺序
    """
    import heapq
    
    degree = [len(G[i]) for i in range(n)]
    order = []
    eliminated = [False] * n
    
    # 使用优先队列按度排序
    heap = [(degree[i], i) for i in range(n)]
    heapq.heapify(heap)
    
    while heap:
        deg, v = heapq.heappop(heap)
        if eliminated[v] or degree[v] != deg:
            continue
            
        order.append(v)
        eliminated[v] = True
        
        # 更新受影响顶点的度
        for u in G[v]:
            if not eliminated[u]:
                # 找到v的邻居中未消去的
                new_neighbors = set()
                for w in G[v]:
                    if not eliminated[w] and w != u:
                        # u和w可能成为新邻居
                        new_neighbors.add(w)
                
                # 更新u的邻接表
                G[u].update(new_neighbors)
                degree[u] = len(G[u])
                heapq.heappush(heap, (degree[u], u))
    
    return order

 

5.3 嵌套分割排序

def nested_dissection(G, n):
    """
    嵌套分割排序：递归地将图分割成两个大致相等的部分
    返回：新的顶点顺序
    """
    if n <= 1:
        return list(range(n))
    
    # 找到图的顶点分割（简化版本）
    separator = find_separator(G, n)
    
    # 将图分成三部分：左、分隔符、右
    left = [v for v in range(n) if v not in separator]
    right = [v for v in range(n) if v in separator]
    
    # 递归排序
    left_order = nested_dissection(restrict_graph(G, left), len(left))
    right_order = nested_dissection(restrict_graph(G, right), len(right))
    
    # 组合顺序：先左，再右，最后分隔符
    return [left[i] for i in left_order] + \
           [right[i] for i in right_order] + \
           list(separator)

 

5.4 并行Cholesky分解

利用消去树进行并行计算：

def parallel_cholesky(A, etree):
    """
    基于消去树的并行Cholesky分解
    etree: 消去树，每个节点记录其父节点和子节点列表
    """
    n = len(A)
    L = np.zeros((n, n))
    
    # 按拓扑顺序处理节点（从叶节点到根）
    def process_node(i):
        # 等待所有子节点完成
        for child in etree.children[i]:
            wait(child)
        
        # 执行当前节点的分解
        # 1. 更新当前行/列（使用子节点的结果）
        for child in etree.children[i]:
            update_from_child(L, i, child)
        
        # 2. 计算当前节点的对角元
        L[i, i] = sqrt(A[i, i] - sum(L[i, :i]**2))
        
        # 3. 计算当前节点所在列的非对角元
        for j in range(i+1, n):
            if A[j, i] != 0 or has_fill(j, i):
                L[j, i] = (A[j, i] - sum(L[j, :i] * L[i, :i])) / L[i, i]
        
        # 标记当前节点完成
        signal_completion(i)
    
    # 并行启动所有叶节点
    for i in range(n):
        if not etree.children[i]:  # 叶节点
            spawn(process_node, i)
    
    return L

 

6. 实际例子

6.1 原始图

考虑矩阵：

对应的图是一条路径 1-2-3-4。

6.2 自然顺序的填充

消元顺序 ：

• 消去1：无填充

• 消去2：2的邻居{1,3}，1和3之间产生填充边(1,3)

• 消去3：3的邻居{2,4}（以及填充的1），产生填充边(2,4)和(1,4)

填充图是完全图K4。

6.3 重排序减少填充

按顺序  ：

• 消去2：邻居{1,3}，产生填充(1,3)

• 消去1：邻居{2,3}（2已消），无新填充

• 消去3：邻居{1,2,4}，1和2已消，只可能产生(1,4)和(2,4)

填充图比完全图少一些边。

7. 算法复杂度

算法	时间复杂度	空间复杂度
符号分解
最小度排序	平均
嵌套分割
并行分解

8. 总结

    图论在Cholesky分解中的应用：

1. 预测填充，即通过符号分解提前知道L的非零模式

2. 减少填充，即通过图重排序算法优化消元顺序

3. 并行计算，利用消去树实现任务并行

4. 内存管理，预分配L的存储空间

5. 数值稳定性，图的信息可用于指导选主元策略

        这些技术是高性能稀疏直接求解器的核心，在有限元分析、电路仿真、优化问题等领域有广泛应用。