MATLAB矩阵：lu、qr、chol 三大分解：解方程、最小二乘、正定矩阵各自的分解武器

一、为什么“分解即降维”：一条贯穿全文的主线

先抛出一个问题：为什么数值线性代数里，几乎所有高效算法都绕不开矩阵分解？笔者的回答是——分解的本质，是把一个“耦合的、难解的原问题”，拆成若干个“解耦的、结构良好的子问题”。这个过程在数学上叫因子化，在工程直觉上就是降维：把一个 n×n 的稠密算子，变成两个三角矩阵（LU）、一个正交矩阵加一个三角矩阵（QR）、或者一个三角矩阵的转置乘积（Cholesky）。

以求解线性方程组 Ax = b 为例。直接求逆 A⁻¹ 再乘 b，理论复杂度是 O(n³)，但数值上极不稳定，且当 b 变化时需要重算。而 LU 分解把 A 拆成 PA = LU 后，求解退化为“前代 + 回代”两个 O(n²) 的三角求解过程。这意味着：一次分解的 O(n³) 成本，可以摊薄到任意多个右端项上。这正是分解“降维”价值的第一个体现——把重复求解的边际成本从立方级压到平方级。

本文评述：很多教材把 LU、QR、Cholesky 并列成三个孤立的知识点，读者学完却不知道何时用哪个。笔者认为，真正有价值的认知框架是：先看矩阵结构（对称？正定？超定？稀疏？），再定分解类型，最后定求解路径。结构决定分解，分解决定算法，算法决定性能。这条链条，就是全文的组织逻辑。

为了把这条主线讲透，本文会反复回到三个维度做对比：数值稳定性（舍入误差如何被放大或抑制）、计算复杂度（时间与空间的量级）、结构利用度（是否吃透了对称性、稀疏性、正定性）。这三个维度，构成了后文所有选型建议的底层依据。

二、LU 分解：通用稠密求解的工业标准

2.1 数学本质：高斯消元法的矩阵语言

LU 分解的数学定义很朴素：对矩阵 A，寻找下三角矩阵 L（对角线为 1）和上三角矩阵 U，使得 A = LU。它的算法内核，其实就是中学就学过的高斯消元——只不过用矩阵乘法把每一步消元“打包”成了初等变换矩阵的乘积。

但朴素 LU 有一个致命弱点：如果主元（pivot）为零或极小，消元过程会产生巨大的舍入误差甚至除零。于是有了带部分选主元的 LU（partial pivoting）：在每一列消元前，把该列绝对值最大的元素换到对角线上。MATLAB 的 lu 函数默认就做部分选主元，返回形式为 PA = LU 或 [L,U,P] = lu(A)。

% 基础 LU 分解与求解
A = [4 -2 1; -2 4 -2; 1 -2 4];
b = [11; -16; 17];

[L, U, P] = lu(A);        % PA = LU
y = L \ (P*b);            % 前代
x = U \ y;                % 回代
disp(x);                  % 验证 Ax = b

% 更推荐：直接反斜杠，MATLAB 自动选最优分解
x2 = A \ b;

这里要特别提醒一个工程习惯：绝大多数情况下，直接写 A \ b 比手动 lu 再前代回代更优。因为 MATLAB 的反斜杠运算符会根据矩阵结构自动选择 LU、QR 或 Cholesky，还会做稀疏性分析。手动分解的价值在于：需要复用分解结果、需要检查条件数、或者需要精细控制数值行为时。

2.2 选主元策略：稳定性与稀疏性的博弈

选主元不是“越稳越好”这么简单。对于稀疏矩阵，部分选主元可能引入大量填充（fill-in），把原本稀疏的 L、U 变得稠密，内存和计算量双双爆炸。这就引出了 LU 分解中最有工程味道的一个话题：如何在数值稳定性和稀疏性之间做权衡。

MATLAB 提供了四种 pivot 选项，这是很多使用者忽略的“隐藏武器”：

pivot 选项	策略	适用场景	风险
'vector'	部分选主元（默认）	通用稠密矩阵	稀疏时填充多
'matrix'	返回置换矩阵 P	需要显式 P 的场合	内存略高
0	不选主元	已知对角占优	可能除零/失稳
'vector' + 重排	先 perm 再分解	稀疏矩阵	需人工调参

对于稀疏矩阵，MATLAB 提供了 decomposition 对象和 symrcm、colamd、dissect 等重排函数。重排的目标很明确：让消元过程中产生的填充最少。这是稀疏直接求解器的核心技巧。

% 稀疏矩阵 LU：重排 + 分解
A = gallery('wathen', 20, 20);   % 生成稀疏对称正定测试矩阵
A = A + speye(size(A));          % 保证正定

p = symrcm(A);                   % 反向 Cuthill-McKee 重排，减小带宽
Aperm = A(p, p);

[L, U, P] = lu(Aperm);           % 对重排后矩阵做 LU
% 求解时记得把 b 也按 p 重排
% x(p) = U \ (L \ (P * b(p)));

笔者认为：稀疏 LU 的工程难度，远大于稠密 LU。稠密 LU 的复杂度是确定的 O(n³)，而稀疏 LU 的填充量高度依赖矩阵的图结构（graph structure）。这也是为什么 MATLAB 的 decomposition 对象值得推荐——它把重排、分解、求解封装成一个可复用的对象，避免每次求解都重新分析结构。

2.3 数值稳定性：增长因子与条件数

衡量 LU 稳定性的经典指标是“增长因子”（growth factor）ρ，定义为消元过程中出现的最大元素与原始最大元素之比。部分选主元能把 ρ 控制在 2ⁿ⁻¹ 的理论上界内，但实际中通常远小于这个界。更实用的判据是条件数 cond(A)：

% 条件数诊断：判断解的可信度
A = hilb(12);                    % 希尔伯特矩阵，著名的病态矩阵
b = ones(12, 1);
x = A \ b;

c = cond(A);                     % 2-范数条件数
r = norm(A*x - b) / norm(b);     % 相对残差
fprintf('cond(A) = %.2e, 相对残差 = %.2e\n', c, r);
% 经验法则：损失约 log10(cond) 位有效数字

这里有一条非常实用的经验法则：如果 cond(A) ≈ 10^k，那么双精度求解大约会损失 k 位有效数字。双精度有约 16 位有效数字，所以当 cond(A) 超过 10¹⁶ 时，解基本不可信。希尔伯特矩阵正是这样一个“照妖镜”，n=12 时条件数已达 10¹⁶ 量级。

三、QR 分解：最小二乘与正交化的中流砥柱

3.1 正交性的魔力：为什么 QR 比正规方程更可靠

最小二乘问题的经典解法是构造正规方程 A^TA x = A^Tb，然后用 Cholesky 或 LU 求解。这条路看起来很美，但有一个致命的数值缺陷：cond(A^TA) = cond(A)²。条件数被平方，意味着有效数字损失翻倍。

QR 分解则完全绕开了这个问题。把 A（m×n，m≥n）分解为 A = QR，其中 Q 是 m×n 正交矩阵（Q^TQ = I），R 是 n×n 上三角矩阵。正交变换不改变 2-范数，因此最小二乘解可以直接由 Rx = Q^Tb 回代得到，条件数保持为 cond(A) 而非其平方。这是 QR 在最小二乘中不可替代的根本原因。

% 最小二乘：正规方程 vs QR，对比精度
m = 200; n = 50;
A = randn(m, n) * 1e-3;          % 制造一定病态性
A(:, end) = A(:, 1) + 1e-10*randn(m,1);  % 制造近共线列
b = randn(m, 1);

% 方法一：正规方程（不推荐）
x_ne = (A'*A) \ (A'*b);

% 方法二：QR（推荐）
[Q, R] = qr(A, 0);               % 经济型 QR
x_qr = R \ (Q'*b);

% 方法三：反斜杠（MATLAB 自动选 QR）
x_bs = A \ b;

fprintf('正规方程残差: %.3e\n', norm(A*x_ne - b));
fprintf('QR 残差:      %.3e\n', norm(A*x_qr - b));
fprintf('反斜杠残差:   %.3e\n', norm(A*x_bs - b));

3.2 三种 QR 算法：Householder、Givens、Gram-Schmidt

QR 分解不是一种算法，而是一族算法。它们在稳定性、并行性、稀疏性上各有取舍：

算法	核心操作	稳定性	最佳场景
Householder	反射变换，一次消一列	优秀	稠密矩阵，MATLAB 默认
Givens	旋转变换，一次消一个元素	优秀	稀疏/带状矩阵，可并行
经典 Gram-Schmidt	逐列正交化	差	教学演示
修正 Gram-Schmidt	投影后立即更新	中等	并行实现，历史遗留代码

MATLAB 的 qr 函数默认使用 Householder 反射，对稀疏矩阵则会自动切换到 Givens 旋转以减少填充。这一点从 qr 对稀疏输入返回的结果结构就能看出。

% 稀疏 QR：观察填充
A = sparse(gallery('poisson', 30));   % 二维泊松矩阵，稀疏
[Q, R] = qr(A);
fprintf('A 非零元: %d\n', nnz(A));
fprintf('R 非零元: %d\n', nnz(R));
fprintf('Q 非零元: %d\n', nnz(Q));
% 稀疏 QR 的 Q 通常比 A 稠密得多，这是填充的代价

本文评述：Gram-Schmidt 家族在教科书里出场率极高，但在生产代码里几乎绝迹。原因很简单：正交性的丢失是不可逆的。经典 Gram-Schmidt 在病态矩阵上会迅速失去列向量的正交性，而 Householder 反射因为每一步都是正交变换，误差不会累积放大。选算法时，稳定性永远优先于实现的简洁性。

3.3 带列主元的 QR：秩亏问题的解药

当 A 秩亏或近秩亏时，普通 QR 得到的 R 对角元会出现极小值，求解不稳定。带列主元的 QR（column pivoting）通过动态交换列，让 R 的对角元按绝对值递减排列，从而暴露矩阵的数值秩。MATLAB 写法是 [Q,R,P] = qr(A)，满足 AP = QR。

% 秩亏最小二乘：列主元 QR
A = [1 2 3; 4 5 6; 7 8 9; 10 11 12];  % 秩为 2 的矩阵
b = [1; 2; 3; 4];

[Q, R, P] = qr(A, 0);
disp(diag(R));                   % 第三个对角元接近 0，暴露秩亏

% 截断求解：只保留前 rank 个分量
tol = max(size(A)) * eps(norm(R));
rank_A = sum(abs(diag(R)) > tol);
x_trunc = R(1:rank_A, 1:rank_A) \ (Q(:, 1:rank_A)' * b);
x_full = zeros(size(A, 2), 1);
x_full(P(1:rank_A)) = x_trunc;
disp(x_full);

四、Cholesky 分解：正定矩阵的专属快车道

4.1 对称正定：一个“特权结构”

如果矩阵 A 对称正定（SPD），那么存在唯一的下三角矩阵 L（对角元为正），使得 A = LL^T。这就是 Cholesky 分解。它相比 LU 有三大优势：

• 计算量减半：约 n³/6 次浮点运算，而 LU 是约 2n³/3，正好是三分之一到一半。
• 存储减半：只需存一个三角矩阵，LU 要存两个。
• 无需选主元：正定性保证了对角元始终为正，消元过程天然稳定。

正因为这些优势，凡是能转化为 SPD 系统的问题，都应该优先走 Cholesky 路线。典型场景包括：有限元刚度矩阵、高斯过程协方差矩阵、岭回归的正则化正规方程、卡尔曼滤波的协方差更新等。

% Cholesky 分解：SPD 系统的标准解法
n = 500;
A = randn(n);
A = A' * A + n * eye(n);        % 构造 SPD 矩阵
b = randn(n, 1);

% 方法一：chol + 两次三角求解
L = chol(A, 'lower');            % A = L*L'
y = L \ b;
x = L' \ y;

% 方法二：反斜杠（MATLAB 自动检测 SPD 并用 Cholesky）
x2 = A \ b;

fprintf('chol 残差: %.3e\n', norm(A*x - b));
fprintf('反斜杠残差: %.3e\n', norm(A*x2 - b));

% 性能对比
f1 = @() chol(A, 'lower') \ b;
f2 = @() A \ b;
t1 = timeit(f1);
t2 = timeit(f2);
fprintf('chol 路径: %.4f ms\n', t1*1000);
fprintf('反斜杠:    %.4f ms\n', t2*1000);

4.2 正定性判据：chol 是“会报错的检测器”

一个常被忽视的实用技巧：用 chol 来检测矩阵是否正定。如果矩阵非正定，chol 会抛出错误或返回非零状态。这比计算所有特征值（O(n³) 且更贵）要高效得多。

% 用 chol 检测正定性
function tf = isPosDef(A)
    try
        chol(A);
        tf = true;
    catch
        tf = false;
    end
end

% 测试
A1 = [2 -1; -1 2];
A2 = [1 2; 2 1];
fprintf('A1 正定? %d\n', isPosDef(A1));   % 1
fprintf('A2 正定? %d\n', isPosDef(A2));   % 0

% 更稳健：带容差的 chol
[R, p] = chol(A2);
if p > 0
    fprintf('非正定，失败于第 %d 个主元\n', p);
end

4.3 半正定与修正 Cholesky

现实中大量矩阵是半正定（PSD）而非严格正定，比如协方差矩阵在样本数小于维度时、图拉普拉斯矩阵在存在零特征值时。此时标准 Cholesky 会失败。工程上有两条主流应对路径：

路径一：对角扰动（jitter）。给对角线加上一个极小正数 ε，把 PSD 变成 SPD。这是高斯过程回归中的标准操作，ε 通常取 1e-6 到 1e-10 之间，需要根据矩阵尺度调整。

路径二：修正 Cholesky（modified Cholesky）。在分解过程中动态检测非正的对角元，将其修正为一个小的正数。Gill-Murray-Wright 算法是经典代表，MATLAB 没有内置，但可以用 ldl 函数配合修正实现。

% 半正定矩阵的 jitter 处理
K = randn(100, 80);
K = K * K';                      % 秩最多 80，半正定
[~, p] = chol(K);
fprintf('原始 chol 失败标志: %d\n', p);   % p > 0，失败

% 加 jitter
jitter = 1e-8 * trace(K) / size(K, 1);
K_reg = K + jitter * eye(size(K));
L = chol(K_reg, 'lower');
fprintf('加 jitter 后成功，jitter = %.2e\n', jitter);

笔者认为：jitter 的选取是一门“玄学”与科学的结合。太小则无法恢复正定性，太大则扭曲原问题。一个可操作的经验是：先取 jitter = 1e-10 × trace(A)/n，若 chol 仍失败，逐步放大 10 倍，直到成功。同时要监控解的变化，若 jitter 超过矩阵尺度的 1e-6，就应警惕结果已偏离原问题。

五、三大分解横向对比与选型决策树

把三大分解放在同一张表里对比，是理解它们分工的最快方式。下表从复杂度、稳定性、结构要求、典型应用四个维度做了汇总（数据来源：Golub & Van Loan《Matrix Computations》第4版，以及 MATLAB R2024b 官方文档）：

维度	LU	QR	Cholesky
矩阵要求	方阵，非奇异	任意 m×n	对称正定
浮点运算量	~2n³/3	~2mn² − 2n³/3	~n³/3
存储	L + U	Q + R	仅 L
稳定性	需选主元	天然稳定	无需选主元
条件数影响	cond(A)	cond(A)	cond(A)
典型应用	通用线性求解	最小二乘、特征值	SPD 系统、优化

基于这张表，笔者提炼出一棵可直接落地的选型决策树：

输入：矩阵 A，右端项 b
│
├─ A 是方阵？
│   ├─ 是 → A 对称？
│   │        ├─ 是 → A 正定？
│   │        │        ├─ 是 → 用 chol（最快最省）
│   │        │        └─ 否 → 用 ldl（对称不定）
│   │        └─ 否 → 用 lu（带选主元）
│   └─ 否（m ≠ n）→ 用 qr（最小二乘）
│
├─ A 稀疏？
│   ├─ 是 → 先重排（symrcm/colamd/dissect）
│   │        再选分解，优先用 decomposition 对象
│   └─ 否 → 直接调用 \, MATLAB 自动选最优
│
└─ 需要复用分解？
    ├─ 是 → 显式分解并保存因子
    └─ 否 → 直接用 A \ b

六、工程实战：从稠密到稀疏、从 CPU 到 GPU

6.1 稠密大规模：分块与内存布局

当矩阵规模达到数千到上万时，性能瓶颈往往不是浮点运算，而是内存带宽。MATLAB 底层调用的是 Intel MKL 或 OpenBLAS，这些库都实现了分块算法（blocked algorithm），把大矩阵切成能装进 L2/L3 缓存的小块，最大化数据复用。用户层面能做的优化有限，但有一点很关键：尽量使用列优先的内存访问模式，因为 MATLAB 是列主序存储。

% 列优先 vs 行优先：性能差异演示
n = 4000;
A = randn(n);

% 好：按列访问
tic;
s1 = 0;
for j = 1:n
    s1 = s1 + sum(A(:, j));
end
t_col = toc;

% 差：按行访问
tic;
s2 = 0;
for i = 1:n
    s2 = s2 + sum(A(i, :));
end
t_row = toc;

fprintf('列访问: %.4f s, 行访问: %.4f s, 比值: %.1fx\n', ...
        t_col, t_row, t_row/t_col);

6.2 稀疏矩阵：重排策略的威力

稀疏矩阵的分解性能，几乎完全由重排策略决定。MATLAB 提供了多种重排算法，它们的适用场景不同：

重排函数	算法	目标	适用结构
symrcm	反向 Cuthill-McKee	减小带宽	带状、对称
symamd	近似最小度	减少填充	对称、不规则
colamd	列近似最小度	减少填充	非对称、超定
dissect	嵌套剖分	分治减少填充	二维/三维网格

% 重排策略对比：以二维泊松矩阵为例
A = gallery('poisson', 50);       % 2500×2500 稀疏矩阵
b = ones(size(A, 1), 1);

% 不重排
tic; x0 = A \ b; t0 = toc;

% symrcm 重排
p1 = symrcm(A);
tic; x1 = A(p1,p1) \ b(p1); t1 = toc;

% symamd 重排
p2 = symamd(A);
tic; x2 = A(p2,p2) \ b(p2); t2 = toc;

% dissect 重排
p3 = dissect(A);
tic; x3 = A(p3,p3) \ b(p3); t3 = toc;

fprintf('无重排:   %.4f s\n', t0);
fprintf('symrcm:   %.4f s\n', t1);
fprintf('symamd:   %.4f s\n', t2);
fprintf('dissect:  %.4f s\n', t3);
% 注：以上为模拟环境下的相对趋势，实际数值随硬件和 MATLAB 版本变化

6.3 GPU 加速：gpuArray 与分解的适配

MATLAB 的 Parallel Computing Toolbox 支持把矩阵搬到 GPU 上做分解。但要注意：不是所有分解都能从 GPU 获得加速。LU 和 Cholesky 在 GPU 上有成熟的 cuSOLVER 实现，加速比可观；而稀疏 QR 的 GPU 支持相对有限，因为稀疏结构导致负载不均衡。

% GPU 加速的 Cholesky 分解
if canUseGPU
    n = 8000;
    A = randn(n, 'single');
    A = A' * A + n * eye(n, 'single');
    b = randn(n, 1, 'single');

    % CPU
    tic; x_cpu = A \ b; t_cpu = toc;

    % GPU
    Ag = gpuArray(A);
    bg = gpuArray(b);
    tic; x_gpu = Ag \ bg; t_gpu = toc;
    x_gpu = gather(x_gpu);

    fprintf('CPU: %.4f s, GPU: %.4f s, 加速比: %.1fx\n', ...
            t_cpu, t_gpu, t_cpu/t_gpu);
    fprintf('结果差异: %.3e\n', norm(x_cpu - x_gpu) / norm(x_cpu));
else
    disp('未检测到可用 GPU');
end

本文评述：GPU 加速的收益高度依赖矩阵规模。经验上，单精度下矩阵维度低于 1000 时，数据传输开销会吃掉全部加速收益。只有当维度达到数千以上、且能保持数据长期驻留 GPU 时，加速才真正划算。盲目上 GPU 是常见的工程误区。

七、前沿预判：混合精度、随机算法与自动微分

7.1 混合精度分解：用低精度换速度

近年来，混合精度（mixed precision）成为数值线性代数的热点。核心思路是：用低精度（如 fp16 或 bf16）做分解的“预测”，再用高精度（fp32/fp64）做迭代精化。这在 GPU 上尤其有效，因为低精度运算的吞吐量可以是高精度的数倍。

2023 年以来，多项研究（如 Higham 团队关于 mixed precision iterative refinement 的系统性工作）表明，在合理设计下，混合精度可以在保持接近双精度精度的同时获得 2-4 倍加速。MATLAB 目前对混合精度的原生支持有限，但在 single 与 double 之间手动实现迭代精化是可行的。

% 混合精度迭代精化：single 分解 + double 修正
n = 3000;
A = randn(n);
A = A' * A + n * eye(n);
b = randn(n, 1);

% 低精度分解
As = single(A);
bs = single(b);
x = double(As \ bs);             % 初始解（低精度）

% 迭代精化
for k = 1:5
    r = b - A * x;               % 高精度残差
    dx = double(As \ single(r)); % 低精度修正
    x = x + dx;
    fprintf('迭代 %d: 残差 = %.3e\n', k, norm(r)/norm(b));
end

7.2 随机数值线性代数：用概率换规模

当矩阵大到无法完整分解时，随机算法提供了另一条路。随机 QR（randomized QR）通过先乘一个随机投影矩阵，把大矩阵压缩成小矩阵，再对小矩阵做标准 QR。这类算法在 2011 年 Halko、Martinsson、Tropp 的综述后被广泛采用，近三年在机器学习和大规模科学计算中应用激增。

MATLAB 没有内置随机 QR，但实现并不复杂。它的价值在于：当矩阵的数值秩远小于维度时，随机 QR 可以把 O(mn²) 的复杂度降到 O(mnk)，k 是目标秩。

% 随机 QR：低秩近似
function [Q, R] = rand_qr(A, k, p)
    % A: m×n, k: 目标秩, p: 过采样参数
    if nargin < 3, p = 10; end
    l = min(k + p, size(A, 2));
    Omega = randn(size(A, 2), l);
    Y = A * Omega;               % 随机投影
    [Q, ~] = qr(Y, 0);           % 正交化
    B = Q' * A;
    [Q2, R] = qr(B, 0);          % 小矩阵 QR
    Q = Q * Q2;
end

% 测试：低秩矩阵
m = 5000; n = 2000; k = 50;
A = randn(m, k) * randn(k, n);   % 数值秩为 k
[Q, R] = rand_qr(A, k);
fprintf('近似误差: %.3e\n', norm(A - Q*R, 'fro') / norm(A, 'fro'));

7.3 自动微分与可微分解

在深度学习与科学计算融合的趋势下，“可微分解”成为一个新兴方向。核心问题是：如何对矩阵分解求导？这在端到端训练包含线性求解层的模型时至关重要。MATLAB 的 Deep Learning Toolbox 和 Automatic Differentiation 支持部分分解的自动求导，但覆盖度仍在完善中。

从数学上看，Cholesky 分解的导数有闭式解（Murray 2016 给出了完整推导），LU 和 QR 的导数则更复杂，需要处理置换和符号问题。笔者认为，这一方向在未来三到五年会持续升温，尤其是当物理信息神经网络（PINN）与经典数值求解器深度耦合时。

八、常见陷阱与调试清单

即便理解了原理，工程实践中仍会踩坑。以下是笔者总结的高频陷阱清单：

1. 用 inv(A)*b 求解：这是最经典的错误。inv 不仅慢（多一次矩阵乘法），而且数值上更不稳定。永远用 A \ b。
2. 忽视条件数：解出来不检查 cond(A)，可能拿到一个完全错误的“解”。建议养成 if cond(A) > 1/eps, warning(...), end 的习惯。
3. 稀疏矩阵用 full 存储：一个 10000×10000 的稀疏矩阵，如果非零元只有 1%，full 存储需要 800MB，sparse 只需 8MB。内存和速度双输。
4. 对非正定矩阵用 chol：chol 会报错，但更隐蔽的是“接近正定”的矩阵，chol 可能成功但结果不可信。建议配合 [~,p] = chol(A) 检查 p。
5. 忘记重排后的索引还原：稀疏矩阵重排后，解向量也需要按逆置换还原，否则结果全错。
6. 在循环里反复分解同一个矩阵：如果右端项变化但矩阵不变，应该分解一次、求解多次。用 decomposition 对象或保存因子。
7. GPU 上做小矩阵分解：数据传输开销大于计算收益，得不偿失。

% 复用分解：decomposition 对象
A = randn(1000);
A = A' * A + 1000 * eye(1000);
dA = decomposition(A, 'chol');    % 一次分解，多次求解

b1 = randn(1000, 1);
b2 = randn(1000, 1);
x1 = dA \ b1;                     % 复用分解
x2 = dA \ b2;
% 比每次 A \ b 快得多

九、结语：把分解当作一种思维方式

回到全文的主线：分解即降维，降维即把难问题转化为易问题的组合。LU 把通用方阵求解降维成两个三角求解；QR 把最小二乘降维成正交投影加回代；Cholesky 把 SPD 系统降维成一次三角分解加两次三角求解。三者共享同一个哲学：识别结构，利用结构，把复杂度压到结构允许的最低水平。

这种思维方式的价值，远超 MATLAB 函数调用本身。当你在设计一个算法、优化一段代码、甚至分析一个业务问题时，都可以问自己：这个问题有没有可利用的结构？能不能分解成几个更简单的子问题？分解后的子问题是不是各自有成熟解法？如果答案是肯定的，那么你已经在用数值线性代数的智慧解决实际问题了。

最后给出一条实用建议：在 MATLAB 里，先信任反斜杠，再在必要时手动分解。反斜杠背后的自动选择逻辑，凝聚了几十年数值分析的研究成果。理解它为什么这样选，比盲目替换成手动分解更有价值。

十、参考文献与拓展资源