25、线性系统求解中的分解方法及代码演变

线性系统求解中的分解方法及代码演变

在求解稠密线性系统的过程中，分解方法是一种基础且有效的手段。下面我们将详细介绍这种方法，并探讨其代码实现随着硬件发展的演变。

分解方法的基本原理

对于线性系统 $Ax = b$，其中 $A$ 是 $n$ 阶非奇异矩阵，分解方法的核心思想是将矩阵 $A$ 分解为更简单矩阵的乘积，从而简化问题的求解。具体来说，可将 $A$ 分解为 $A = LU$ 的形式，其中 $L$ 是下三角矩阵（对角线上方元素均为 0，且对角线元素为 1），$U$ 是上三角矩阵（对角线下方元素均为 0）。

在分解过程中，会使用矩阵 $A$ 的对角元素（称为主元）来除对角线下方的元素。若矩阵 $A$ 存在零主元，分解过程会因除零错误而中断；主元值过小也会过度放大计算过程中的数值误差。为保证数值稳定性，需要交换矩阵的行，确保主元的绝对值尽可能大。引入行置换矩阵 $P$ 后，分解形式变为 $P^TA = LU$，此时方程组的解可表示为 $x = A^{-1}Pb$，对应的求解算法步骤如下：
1. 对矩阵 $A$ 进行分解。
2. 求解系统 $Ly = Pb$。
3. 求解系统 $Ux = y$。

这种通过分解进行矩阵计算的方法具有显著优势。它将计算过程分为两个阶段：先进行分解计算，再利用分解结果求解当前问题。若在计算过程中需要求解不同的右侧向量，矩阵 $A$ 只需分解一次，后续可重复使用分解结果。其中，步骤 1 对矩阵 $A$ 进行分解需要 $O(n^3)$ 次运算，而步骤 2 和 3 的求解过程仅需 $O(n^2)$ 次运算。此外，$L$ 和 $U$ 因子无需额外存储空间，可直接占用矩阵 $A$ 原本的空间。