分块矩阵的Tikhonov正则化在病态线性方程组求解中的应用

题目描述
考虑求解病态线性方程组 Ax = b，其中 A 是大型稀疏或稠密的分块矩阵（如来自图像重建或逆问题的离散化），b 是已知向量。由于问题本身的不适定性或矩阵 A 的条件数很大，直接求解可能对数据误差极其敏感，导致解不稳定。Tikhonov 正则化通过引入正则化项来稳定求解，将原问题转化为最小化问题：min{ ||Ax - b||² + λ²||Lx||² }，其中 λ 是正则化参数，控制拟合程度与解平滑性的平衡，L 通常是单位矩阵或离散微分算子（也可能是分块结构）。当 A 和 L 为分块矩阵时，如何高效求解对应的正则化系统是本问题的核心。

解题过程

1. 问题转化：正则化正规方程组
   原始最小化问题 min{ ||Ax - b||² + λ²||Lx||² } 的解等价于求解线性系统：
   (AᵀA + λ²LᵀL) x = Aᵀb。
   这个系统称为正则化正规方程组。由于 AᵀA + λ²LᵀL 是对称正定（或半正定）的，理论上可用 Cholesky 分解或共轭梯度法求解。但当 A 和 L 是大型分块矩阵时，直接形成 AᵀA 或 LᵀL 可能计算昂贵且破坏稀疏性，需寻求更高效方法。

2. 增广系统方法
   为避免显式形成 AᵀA，将问题改写为增广形式：

\[\begin{bmatrix} A \\ \lambda L \end{bmatrix} x \approx \begin{bmatrix} b \\ 0 \end{bmatrix} \]

这等价于求解最小二乘问题：min ||Bx - c||²，其中 B = [A; λL]（纵向分块），c = [b; 0]。此时，可直接对增广矩阵 B 应用 QR 分解或迭代法（如 LSQR）。若 A 和 L 有分块结构（如分块对角或带状），B 会继承类似结构，利于分块运算。

3. 分块矩阵情况下的迭代求解
   当 A 和 L 为分块矩阵（例如来自区域分解离散化）时，常用迭代法如共轭梯度法（CG）求解正规方程。关键步骤是设计高效的矩阵-向量乘法：

   • 对于向量 v，计算 (AᵀA + λ²LᵀL)v 时，避免显式构造 AᵀA，而是依次计算：
     • u₁ = Av（利用 A 的分块结构并行计算），
     • u₂ = Aᵀu₁（分块转置乘法），
     • 类似处理 LᵀLv，
       最后组合 u₂ + λ²(LᵀLv)。
       此方法保持矩阵的稀疏性和分块特性，减少存储和计算量。

4. 正则化参数 λ 的选择
   Tikhonov 正则化的有效性依赖于 λ。常用选择方法包括：

   • L-曲线准则：绘制 log||Ax - b|| 与 log||Lx|| 的曲线，选择拐点对应的 λ。
   • 广义交叉验证（GCV）：最小化 G(λ) = ||Ax - b||² / [trace(I - A(AᵀA + λ²LᵀL)⁻¹Aᵀ)]²。
     在分块矩阵场景中，参数选择可能需结合分块结构简化 trace 计算（如使用随机迹估计）。

5. 预处理技术加速收敛
   病态系统即使正则化后仍可能迭代缓慢。针对分块结构，可设计分块预处理子：

   • 若 L 是单位矩阵，预处理子可取 (AᵀA + λ²I) 的近似，例如用分块对角部分或不完全 Cholesky 分解。
   • 对于一般 L，可基于分块 Schur 补或区域分解构造预处理子，降低条件数。

通过结合分块矩阵运算和迭代法，Tikhonov 正则化能稳定求解病态问题，同时利用结构提升效率。