跳到主要内容

高斯消去

一、直观理解：把方程组“三角化”再回代​

高斯消去法（Gaussian elimination）通过初等行变换，把线性方程组 Ax=bAx=b 的增广矩阵化为上三角形式，再从最后一个未知数开始回代求解。它对应我们手算消元的系统化版本，也是 LU 分解与众多直接法的基础。

直观上，每一步用某一行消去其下方某列的系数，像层层削去变量依赖，直到只剩一个单变量方程。

二、适用场景与优缺点​

类型说明
适用场景中小规模稠密线性系统；教学与原理实现；矩阵结构无特殊可利用时的基线直接法。
核心优点有限步得到精确解（在精确算术下）；复杂度明确 O(n3)O(n^3)；可同时处理多个右端项（先分解后回代）。
局限与缺点无主元时可能数值不稳或中断；稠密填充；大型稀疏问题应改用稀疏直接法或迭代法。

三、核心原理：消元 + 回代​

对 k=1,…,n−1k=1,\ldots,n-1，若主元 akk≠0a_{kk}\ne0，取乘数

ℓik=aikakk,i=k+1,…,n\ell_{ik}=\frac{a_{ik}}{a_{kk}},\quad i=k+1,\ldots,n

第 ii 行减去 ℓik\ell_{ik} 倍第 kk 行，同时对 bb 做相同操作。完成后 AA 成为上三角 UU，再解 Ux=b′Ux=b'。

部分主元（Partial Pivoting）​

在第 kk 列的 k…nk\ldots n 行中选取绝对值最大元作为主元并交换行。这能避免除以过小主元，是实际计算中的标准稳定化手段。MATLAB 的 lu 与 \ 都隐含主元策略。

复杂度与误差​

  • 消元约 23n3\frac{2}{3}n^3 次浮点运算；
  • 误差受 κ(A)=∥A∥∥A−1∥\kappa(A)=\|A\|\|A^{-1}\| 控制：残差小不等于解误差小。

四、实践要点​

  1. 永远配合主元选择；
  2. 用残差 r=b−Axr=b-Ax 与相对残差检查；
  3. 多右端项时不要反复全消元，应保存分解；
  4. 检查条件数 cond(A)。

五、MATLAB 示例​

%% 高斯消去（带部分主元）求解 Ax=b
A = [2 1 -1;
     -3 -1 2;
     -2 1 2];
b = [8; -11; -3];

n = size(A,1);
Ab = [A, b];

for k = 1:n-1
    % 部分主元
    [~, p] = max(abs(Ab(k:n, k)));
    p = p + k - 1;
    Ab([k p], :) = Ab([p k], :);
    if abs(Ab(k,k)) < 1e-14
        error("主元过小，矩阵可能奇异");
    end
    for i = k+1:n
        m = Ab(i,k) / Ab(k,k);
        Ab(i,k:end) = Ab(i,k:end) - m * Ab(k,k:end);
    end
end

% 回代
x = zeros(n,1);
for i = n:-1:1
    x(i) = (Ab(i,end) - Ab(i,i+1:n)*x(i+1:n)) / Ab(i,i);
end

xMatlab = A \ b;
fprintf("手写消元解：\n"); disp(x');
fprintf("MATLAB 反斜杠：\n"); disp(xMatlab');
fprintf("残差范数 = %.3e\n", norm(A*x - b));

病态示例：

%% Hilbert 矩阵：消元能算，但解不可靠
n = 10;
H = hilb(n);
xTrue = ones(n,1);
b = H * xTrue;
x = H \ b;
fprintf("cond(H) = %.3e\n", cond(H));
fprintf("相对误差 = %.3e\n", norm(x-xTrue)/norm(xTrue));
fprintf("相对残差 = %.3e\n", norm(H*x-b)/norm(b));

六、如何解释结果​

报告矩阵规模、是否主元、cond(A)、解向量、残差范数。若条件数很大，应说明解的有效数字可能很少，并考虑重整标度、改模型或正则化。

七、常见误区​

1. 无主元消元​

主元为 0 会中断；主元很小会放大舍入误差。

2. 残差很小就认为解很准​

病态时残差小仍可能解误差大。

3. 把 inv(A)*b 当标准解法​

更慢更不稳，应使用 A\b。

4. 对大型稀疏矩阵用稠密消元​

填充与内存会爆炸。

附：建模流程中的落地清单（高斯消去）​

在竞赛或课程论文中，仅写出定义式往往不够。建议按下列顺序组织材料，既能保证可复现，也便于评委快速抓住贡献。

1. 问题对齐​

先用两三句话说明：题目中的哪个子问题对应本方法，输入是什么，输出要支持什么决策。避免“因为常见所以使用”。若存在更简单的替代方法，应说明为何不够用。

2. 符号与单位​

列出核心符号表，强调单位一致。许多错误并非算法失败，而是日/年、米/千米、概率/百分比混用。对无量纲参数，说明取值来源（估计、文献、情景设定）。

3. 关键公式回顾​

核心是把 Ax=bAx=b 经行变换化为上三角后回代；稳定实现必须配合部分主元。条件数 κ(A)\kappa(A) 决定解对数据扰动的放大倍数。

4. 计算与实现要点​

  • 固定随机种子（若涉及随机）；
  • 记录中间量以便手算核对小例子；
  • 用极限情形或解析特例做冒烟测试；
  • 保存主要输入输出到表格，便于敏感性分析复用。

5. 补充 MATLAB 片段​

% 残差与相对指标
A = gallery('lehmer',6); b = sum(A,2); x = A\b;
r = b-A*x; fprintf('relres=%.3e cond=%.3e\n', norm(r)/norm(b), cond(A));

6. 结果报告模板​

建议至少报告：矩阵规模、主元策略、条件数、解向量、相对残差，以及是否与 A\b 一致。

7. 何时不该使用​

矩阵严重病态却仍把数值解当精确真理；或对超大规模稀疏系统仍用稠密消元。

完成以上清单后，再进入“小结”中的一句话金句，通常论文的方法段与验证段就会扎实许多。

八、小结​

高斯消去通过三角化与回代在有限步内求解线性方程组，实际计算必须配合主元选择，并用条件数与残差共同评价结果。它是理解 LU、主元与直接法稳定性的入口。