随机化Chebyshev迭代法求解大型稀疏线性方程组
字数 6428 2025-12-24 11:45:54

随机化Chebyshev迭代法求解大型稀疏线性方程组

问题描述

我们考虑求解大型稀疏线性方程组：

\[A\mathbf{x} = \mathbf{b}, \]

其中 \(A \in \mathbb{R}^{n \times n}\) 是一个大型稀疏矩阵（即绝大多数元素为零），\(\mathbf{b} \in \mathbb{R}^n\) 是已知的右端向量，\(\mathbf{x} \in \mathbb{R}^n\) 是未知解向量。

在科学计算和工程应用中（如有限元法、计算流体力学等），此类方程组经常出现。当矩阵 \(A\) 规模极大（例如 \(n > 10^6\)）时，直接法（如LU分解）由于内存和计算成本过高而变得不可行。迭代法，特别是Krylov子空间方法（如共轭梯度法CG、广义最小残差法GMRES）成为主流选择。然而，对于高度病态（即条件数 \(\kappa(A)\) 非常大）或特征值分布复杂的矩阵，经典Krylov方法的收敛可能非常缓慢。

随机化Chebyshev迭代法是一种基于多项式加速和随机采样的迭代方法。其核心思想是：通过随机化技术，构造一个Chebyshev多项式，该多项式在 \(A\) 的特征值区间上能够有效压缩残差，从而加速收敛。该方法尤其适用于对称正定矩阵，但也可通过预处理或变换扩展到一般矩阵。与确定性Chebyshev迭代相比，随机化版本通过引入随机性来避免精确计算特征值边界，从而降低计算开销。

核心思路

假设 \(A\) 是对称正定的（若不对称，可考虑正规化或预处理后近似对称正定）。经典Chebyshev迭代需要已知 \(A\) 的特征值区间 \([ \lambda_{\min}, \lambda_{\max} ]\)，并利用Chebyshev多项式在该区间上的极小极大性质来加速收敛。但精确计算 \(\lambda_{\min}\) 和 \(\lambda_{\max}\) 成本高昂。

随机化Chebyshev迭代法的创新点在于：

  1. 使用随机化方法（如幂迭代或Lanczos的随机变体）快速估计特征值边界。
  2. 基于估计的边界构造Chebyshev多项式，进行迭代求解。
  3. 在迭代过程中，可能结合随机采样技术来降低每次迭代的计算量（例如，随机选取部分行或列进行计算）。

下面，我们逐步详细讲解算法的推导与实现。

第一步：回顾经典Chebyshev迭代法（确定性版本）

对于对称正定矩阵 \(A\)，其条件数 \(\kappa = \frac{\lambda_{\max}}{\lambda_{\min}}\)。我们希望求解 \(A\mathbf{x} = \mathbf{b}\)。

Chebyshev迭代的格式为：

\[\mathbf{x}_{k+1} = \mathbf{x}_k + \alpha_k (\mathbf{b} - A\mathbf{x}_k), \]

但这里 \(\alpha_k\) 不是常数，而是根据Chebyshev多项式的递推关系设计。

更常见的形式是三项递推：

\[\mathbf{x}_{k+1} = \omega_{k+1} \left( \mathbf{b} - A\mathbf{x}_k + \beta_k \mathbf{x}_k \right) + (1 - \omega_{k+1}) \mathbf{x}_{k-1}, \]

其中参数 \(\omega_{k+1}, \beta_k\) 依赖于特征值区间 \([a, b]\)（这里 \(a = \lambda_{\min}, b = \lambda_{\max}\)）。

具体地，设 \(\gamma = \frac{b - a}{2}, \delta = \frac{b + a}{2}\)，则Chebyshev多项式 \(T_k(t)\) 在区间 \([a, b]\) 上缩放后，其递推为：

\[\mathbf{x}_{k+1} = \frac{2}{\gamma} \left( \frac{A - \delta I}{\gamma} \right) \mathbf{x}_k - \mathbf{x}_{k-1}. \]

但这需要归一化。实用中，采用以下参数：

\[\omega_1 = \frac{2}{b}, \quad \omega_{k+1} = \frac{4}{b - a} \cdot \frac{1}{2 - \frac{b+a}{2} \omega_k}, \quad \beta_k = \frac{b - a}{4} \omega_k. \]

迭代格式为：

\[\mathbf{x}_{k+1} = \mathbf{x}_k + \omega_{k+1} (\mathbf{b} - A\mathbf{x}_k) + \beta_k (\mathbf{x}_k - \mathbf{x}_{k-1}). \]

该方法的收敛速度依赖于条件数 \(\kappa\)：残差 \(\|\mathbf{b} - A\mathbf{x}_k\|\) 以 \(O\left( \left( \frac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1} \right)^k \right)\) 衰减。

问题：需要已知 \(a, b\)（即 \(\lambda_{\min}, \lambda_{\max}\)）。对于大型稀疏矩阵，精确计算它们代价很大。

第二步：随机化特征值边界估计

我们不需要高精度边界，只需足够好的估计来构造多项式。随机化方法可以快速给出估计。

方法1：随机幂迭代（Randomized Power Method）

  • 生成一个随机向量 \(\mathbf{v}_0 \in \mathbb{R}^n\)（通常服从标准正态分布）。
  • 迭代：\(\mathbf{v}_{k+1} = \frac{A \mathbf{v}_k}{\|A \mathbf{v}_k\|}\)。
  • 经过 \(m\) 步后，Rayleigh商 \(\mathbf{v}_m^\top A \mathbf{v}_m\) 给出 \(\lambda_{\max}\) 的近似。
  • 为了估计 \(\lambda_{\min}\)，可以对 \(A^{-1}\) 做幂迭代（但需求解线性系统，不现实），或对 \(A\) 做幂迭代求最大特征值，然后利用矩阵平移：例如，计算 \(B = \sigma I - A\) 的最大特征值，其中 \(\sigma\) 是 \(A\) 的某个上界估计（如 \(\|A\|_\infty\)），则 \(\lambda_{\min}(A) \approx \sigma - \lambda_{\max}(B)\)。

方法2：随机Lanczos（Randomized Lanczos）

  • 更高效，可同时估计最大和最小特征值。
  • 从随机向量 \(\mathbf{v}_1\) 开始，运行 \(m\) 步Lanczos过程，得到三对角矩阵 \(T_m\)。
  • 计算 \(T_m\) 的特征值，它们近似 \(A\) 的极端特征值。
  • 取 \(T_m\) 的最大特征值作为 \(\tilde{\lambda}_{\max}\)，最小特征值作为 \(\tilde{\lambda}_{\min}\)。
  • 通常 \(m \approx 20 \sim 50\) 已足够给出较好的估计。

关键：随机化方法以高概率给出正确估计，且计算成本 \(O(m \cdot \text{nnz}(A))\)（其中 \(\text{nnz}(A)\) 是 \(A\) 的非零元数），远低于精确计算。

第三步：随机化Chebyshev迭代算法

假设通过随机化方法得到特征值边界估计 \(\tilde{a} \approx \lambda_{\min}\)，\(\tilde{b} \approx \lambda_{\max}\)。我们设置一个安全边际（safety margin）以防止估计偏差导致算法不稳定：例如，取 \(a = \tilde{a} - \epsilon\)，\(b = \tilde{b} + \epsilon\)，其中 \(\epsilon\) 是一个小的正数（如 \(0.1 \times (\tilde{b} - \tilde{a})\)）。

然后，使用经典Chebyshev迭代公式，但基于 \([a, b]\) 而非真实区间。迭代步骤如下：

  1. 初始化：

    • 选择初始猜测 \(\mathbf{x}_0\)（通常设为零向量或随机向量）。
    • 计算初始残差 \(\mathbf{r}_0 = \mathbf{b} - A\mathbf{x}_0\)。
    • 设 \(\mathbf{x}_1 = \mathbf{x}_0 + \omega_1 \mathbf{r}_0\)，其中 \(\omega_1 = \frac{2}{b}\)（这是Chebyshev迭代的第一步特殊形式）。
    • 计算残差 \(\mathbf{r}_1 = \mathbf{b} - A\mathbf{x}_1\)。

  2. 迭代（对于 \(k = 1, 2, \dots, K_{\max}\)）：

    • 计算参数：

\[ \omega_{k+1} = \frac{4}{b - a} \cdot \frac{1}{2 - \frac{b+a}{2} \omega_k}, \quad \beta_k = \frac{b - a}{4} \omega_k. \]

  • 更新解：

\[ \mathbf{x}_{k+1} = \mathbf{x}_k + \omega_{k+1} \mathbf{r}_k + \beta_k (\mathbf{x}_k - \mathbf{x}_{k-1}). \]

  • 计算新残差：\(\mathbf{r}_{k+1} = \mathbf{b} - A\mathbf{x}_{k+1}\)。
  • 检查收敛条件：若 \(\|\mathbf{r}_{k+1}\| < \text{tol} \cdot \|\mathbf{b}\|\)，则停止。

注意：由于边界是估计的，实际收敛可能略慢于使用真实边界的Chebyshev迭代，但避免了计算真实边界的昂贵成本。

第四步：随机化采样加速（可选扩展）

对于超大规模问题，即使每次迭代的矩阵-向量乘法 \(A\mathbf{x}\) 也可能很昂贵。我们可以引入随机采样来进一步降低每次迭代的成本。

基本思想：在每次迭代中，不是完整计算 \(A\mathbf{x}_k\)，而是通过随机选择矩阵的行或列的子集来近似。

例如，设 \(S \in \mathbb{R}^{s \times n}\) 是一个随机采样矩阵（如高斯随机矩阵、随机投影或稀疏嵌入矩阵），其中 \(s \ll n\)。则我们用 \(S A S^\top\) 的近似来替代 \(A\) 在局部计算。更具体地，可以将迭代格式改写为：

\[\mathbf{x}_{k+1} = \mathbf{x}_k + \omega_{k+1} (S\mathbf{b} - (S A) \mathbf{x}_k) + \beta_k (\mathbf{x}_k - \mathbf{x}_{k-1}), \]

但需要注意，这种采样会引入额外误差，因此需要权衡计算加速与收敛速度下降。

实践中，更常见的是将随机采样用于预处理或近似矩阵-向量积，而不是直接替代 \(A\)。例如，使用随机化方法构造一个预处理子 \(M \approx A^{-1}\)，然后求解 \(M A \mathbf{x} = M \mathbf{b}\)，此时矩阵 \(M A\) 的条件数改善，Chebyshev迭代收敛更快。

第五步：算法总结与复杂度分析

算法步骤：

  1. 随机化边界估计：

    • 使用随机Lanczos（推荐）或随机幂迭代，运行 \(m\) 步，得到特征值边界估计 \(\tilde{a}, \tilde{b}\)。
    • 加入安全边际：\(a = \tilde{a} - \epsilon, \quad b = \tilde{b} + \epsilon\)。

  2. Chebyshev迭代求解：

    • 初始化 \(\mathbf{x}_0, \mathbf{r}_0\)。
    • 计算第一步：\(\mathbf{x}_1 = \mathbf{x}_0 + \frac{2}{b} \mathbf{r}_0\)。
    • 对于 \(k = 1, 2, \dots\) 直到收敛：
      • 计算 \(\omega_{k+1}, \beta_k\) 如上。
      • 更新 \(\mathbf{x}_{k+1}\) 和残差 \(\mathbf{r}_{k+1}\)。
      • 检查收敛。

计算复杂度：

  • 边界估计：\(O(m \cdot \text{nnz}(A))\) 次浮点运算。
  • 每轮迭代：一次矩阵-向量乘法（\(O(\text{nnz}(A))\)）和几次向量操作（\(O(n)\)）。
  • 总成本：边界估计成本 + 迭代次数 × \(O(\text{nnz}(A))\)。

优势：

  • 避免了精确特征值计算。
  • 对于对称正定矩阵，收敛速度有保证（依赖于估计边界的质量）。
  • 易于并行化（矩阵-向量乘法和向量操作均可并行）。

局限性：

  • 主要适用于对称正定矩阵。对于非对称矩阵，需要转化为正规方程 \(A^\top A \mathbf{x} = A^\top \mathbf{b}\) 或使用其他多项式加速方法。
  • 估计边界不准确时，可能收敛变慢甚至不稳定（安全边际可缓解）。

第六步：数值示例（概念性说明）

假设 \(A\) 是 \(1000 \times 1000\) 的对称正定稀疏矩阵，条件数 \(\kappa \approx 10^3\)。我们想解 \(A\mathbf{x} = \mathbf{b}\)。

  1. 运行随机Lanczos 30步，得到 \(\tilde{\lambda}_{\min} \approx 0.95\)，\(\tilde{\lambda}_{\max} \approx 950\)。
  2. 设置安全边际 \(\epsilon = 5\)，得到区间 \([a, b] = [0.9, 955]\)。
  3. 运行Chebyshev迭代，参数 \(\omega_k, \beta_k\) 根据上述区间计算。
  4. 迭代约 \(k = 50\) 步后，残差下降至 \(10^{-6}\)（相对于初始残差）。

若使用经典CG方法，收敛步数约为 \(O(\sqrt{\kappa}) \approx 32\) 步。随机化Chebyshev可能略多几步（因边界估计误差），但避免了计算精确特征值，总体时间可能更优。

总结

随机化Chebyshev迭代法将经典多项式加速与现代随机化技术结合，为大型稀疏对称正定线性方程组提供了一种高效求解器。其核心贡献在于通过随机化快速估计特征值边界，从而避免了传统Chebyshev迭代中对精确边界的需求。该方法在机器学习和科学计算中具有应用潜力，特别是当矩阵条件数较大且直接求解特征值不可行时。

值得注意的是，该算法可进一步与随机采样、预处理技术结合，以适应更广泛的矩阵类型和计算环境。

随机化Chebyshev迭代法求解大型稀疏线性方程组 问题描述 我们考虑求解大型稀疏线性方程组： \[ A\mathbf{x} = \mathbf{b}, \] 其中 \( A \in \mathbb{R}^{n \times n} \) 是一个大型稀疏矩阵（即绝大多数元素为零），\(\mathbf{b} \in \mathbb{R}^n\) 是已知的右端向量，\(\mathbf{x} \in \mathbb{R}^n\) 是未知解向量。 在科学计算和工程应用中（如有限元法、计算流体力学等），此类方程组经常出现。当矩阵 \(A\) 规模极大（例如 \(n > 10^6\)）时，直接法（如LU分解）由于内存和计算成本过高而变得不可行。迭代法，特别是 Krylov子空间方法 （如共轭梯度法CG、广义最小残差法GMRES）成为主流选择。然而，对于 高度病态 （即条件数 \(\kappa(A)\) 非常大）或 特征值分布复杂 的矩阵，经典Krylov方法的收敛可能非常缓慢。 随机化Chebyshev迭代法 是一种基于 多项式加速 和 随机采样 的迭代方法。其核心思想是：通过随机化技术，构造一个 Chebyshev多项式 ，该多项式在 \(A\) 的特征值区间上能够有效压缩残差，从而加速收敛。该方法尤其适用于 对称正定 矩阵，但也可通过预处理或变换扩展到一般矩阵。与确定性Chebyshev迭代相比，随机化版本通过引入随机性来 避免精确计算特征值边界 ，从而降低计算开销。 核心思路 假设 \(A\) 是对称正定的（若不对称，可考虑正规化或预处理后近似对称正定）。经典Chebyshev迭代需要已知 \(A\) 的特征值区间 \([ \lambda_ {\min}, \lambda_ {\max} ]\)，并利用Chebyshev多项式在该区间上的 极小极大性质 来加速收敛。但精确计算 \(\lambda_ {\min}\) 和 \(\lambda_ {\max}\) 成本高昂。 随机化Chebyshev迭代法的创新点在于： 使用 随机化方法 （如幂迭代或Lanczos的随机变体）快速估计特征值边界。 基于估计的边界构造Chebyshev多项式，进行迭代求解。 在迭代过程中，可能结合 随机采样技术 来降低每次迭代的计算量（例如，随机选取部分行或列进行计算）。 下面，我们逐步详细讲解算法的推导与实现。 第一步：回顾经典Chebyshev迭代法（确定性版本） 对于对称正定矩阵 \(A\)，其条件数 \(\kappa = \frac{\lambda_ {\max}}{\lambda_ {\min}}\)。我们希望求解 \(A\mathbf{x} = \mathbf{b}\)。 Chebyshev迭代的格式为： \[ \mathbf{x}_ {k+1} = \mathbf{x}_ k + \alpha_ k (\mathbf{b} - A\mathbf{x}_ k), \] 但这里 \(\alpha_ k\) 不是常数，而是根据Chebyshev多项式的递推关系设计。 更常见的形式是三项递推： \[ \mathbf{x} {k+1} = \omega {k+1} \left( \mathbf{b} - A\mathbf{x} k + \beta_ k \mathbf{x} k \right) + (1 - \omega {k+1}) \mathbf{x} {k-1}, \] 其中参数 \(\omega_ {k+1}, \beta_ k\) 依赖于特征值区间 \([ a, b]\)（这里 \(a = \lambda_ {\min}, b = \lambda_ {\max}\)）。 具体地，设 \(\gamma = \frac{b - a}{2}, \delta = \frac{b + a}{2}\)，则Chebyshev多项式 \(T_ k(t)\) 在区间 \([ a, b ]\) 上缩放后，其递推为： \[ \mathbf{x} {k+1} = \frac{2}{\gamma} \left( \frac{A - \delta I}{\gamma} \right) \mathbf{x} k - \mathbf{x} {k-1}. \] 但这需要归一化。实用中，采用以下参数： \[ \omega_ 1 = \frac{2}{b}, \quad \omega {k+1} = \frac{4}{b - a} \cdot \frac{1}{2 - \frac{b+a}{2} \omega_ k}, \quad \beta_ k = \frac{b - a}{4} \omega_ k. \] 迭代格式为： \[ \mathbf{x}_ {k+1} = \mathbf{x} k + \omega {k+1} (\mathbf{b} - A\mathbf{x}_ k) + \beta_ k (\mathbf{x} k - \mathbf{x} {k-1}). \] 该方法的收敛速度依赖于条件数 \(\kappa\)：残差 \(\|\mathbf{b} - A\mathbf{x}_ k\|\) 以 \(O\left( \left( \frac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1} \right)^k \right)\) 衰减。 问题 ：需要已知 \(a, b\)（即 \(\lambda_ {\min}, \lambda_ {\max}\)）。对于大型稀疏矩阵，精确计算它们代价很大。 第二步：随机化特征值边界估计 我们不需要高精度边界，只需足够好的估计来构造多项式。随机化方法可以快速给出估计。 方法1：随机幂迭代（Randomized Power Method） 生成一个随机向量 \(\mathbf{v}_ 0 \in \mathbb{R}^n\)（通常服从标准正态分布）。 迭代：\(\mathbf{v}_ {k+1} = \frac{A \mathbf{v}_ k}{\|A \mathbf{v}_ k\|}\)。 经过 \(m\) 步后，Rayleigh商 \(\mathbf{v}_ m^\top A \mathbf{v} m\) 给出 \(\lambda {\max}\) 的近似。 为了估计 \(\lambda_ {\min}\)，可以对 \(A^{-1}\) 做幂迭代（但需求解线性系统，不现实），或对 \(A\) 做幂迭代求最大特征值，然后利用矩阵平移：例如，计算 \(B = \sigma I - A\) 的最大特征值，其中 \(\sigma\) 是 \(A\) 的某个上界估计（如 \(\|A\| \infty\)），则 \(\lambda {\min}(A) \approx \sigma - \lambda_ {\max}(B)\)。 方法2：随机Lanczos（Randomized Lanczos） 更高效，可同时估计最大和最小特征值。 从随机向量 \(\mathbf{v}_ 1\) 开始，运行 \(m\) 步Lanczos过程，得到三对角矩阵 \(T_ m\)。 计算 \(T_ m\) 的特征值，它们近似 \(A\) 的极端特征值。 取 \(T_ m\) 的最大特征值作为 \(\tilde{\lambda} {\max}\)，最小特征值作为 \(\tilde{\lambda} {\min}\)。 通常 \(m \approx 20 \sim 50\) 已足够给出较好的估计。 关键 ：随机化方法以 高概率 给出正确估计，且计算成本 \(O(m \cdot \text{nnz}(A))\)（其中 \(\text{nnz}(A)\) 是 \(A\) 的非零元数），远低于精确计算。 第三步：随机化Chebyshev迭代算法 假设通过随机化方法得到特征值边界估计 \(\tilde{a} \approx \lambda_ {\min}\)，\(\tilde{b} \approx \lambda_ {\max}\)。我们设置一个 安全边际 （safety margin）以防止估计偏差导致算法不稳定：例如，取 \(a = \tilde{a} - \epsilon\)，\(b = \tilde{b} + \epsilon\)，其中 \(\epsilon\) 是一个小的正数（如 \(0.1 \times (\tilde{b} - \tilde{a})\)）。 然后，使用经典Chebyshev迭代公式，但基于 \([ a, b ]\) 而非真实区间。迭代步骤如下： 初始化 ： 选择初始猜测 \(\mathbf{x}_ 0\)（通常设为零向量或随机向量）。 计算初始残差 \(\mathbf{r}_ 0 = \mathbf{b} - A\mathbf{x}_ 0\)。 设 \(\mathbf{x}_ 1 = \mathbf{x}_ 0 + \omega_ 1 \mathbf{r}_ 0\)，其中 \(\omega_ 1 = \frac{2}{b}\)（这是Chebyshev迭代的第一步特殊形式）。 计算残差 \(\mathbf{r}_ 1 = \mathbf{b} - A\mathbf{x}_ 1\)。 迭代（对于 \(k = 1, 2, \dots, K_ {\max}\)） ： 计算参数： \[ \omega_ {k+1} = \frac{4}{b - a} \cdot \frac{1}{2 - \frac{b+a}{2} \omega_ k}, \quad \beta_ k = \frac{b - a}{4} \omega_ k. \] 更新解： \[ \mathbf{x}_ {k+1} = \mathbf{x} k + \omega {k+1} \mathbf{r}_ k + \beta_ k (\mathbf{x} k - \mathbf{x} {k-1}). \] 计算新残差：\(\mathbf{r} {k+1} = \mathbf{b} - A\mathbf{x} {k+1}\)。 检查收敛条件：若 \(\|\mathbf{r}_ {k+1}\| < \text{tol} \cdot \|\mathbf{b}\|\)，则停止。 注意 ：由于边界是估计的，实际收敛可能略慢于使用真实边界的Chebyshev迭代，但避免了计算真实边界的昂贵成本。 第四步：随机化采样加速（可选扩展） 对于超大规模问题，即使每次迭代的矩阵-向量乘法 \(A\mathbf{x}\) 也可能很昂贵。我们可以引入 随机采样 来进一步降低每次迭代的成本。 基本思想：在每次迭代中，不是完整计算 \(A\mathbf{x}_ k\)，而是通过随机选择矩阵的行或列的子集来近似。 例如，设 \(S \in \mathbb{R}^{s \times n}\) 是一个随机采样矩阵（如高斯随机矩阵、随机投影或稀疏嵌入矩阵），其中 \(s \ll n\)。则我们用 \(S A S^\top\) 的近似来替代 \(A\) 在局部计算。更具体地，可以将迭代格式改写为： \[ \mathbf{x}_ {k+1} = \mathbf{x} k + \omega {k+1} (S\mathbf{b} - (S A) \mathbf{x}_ k) + \beta_ k (\mathbf{x} k - \mathbf{x} {k-1}), \] 但需要注意，这种采样会引入额外误差，因此需要权衡计算加速与收敛速度下降。 实践中，更常见的是将随机采样用于 预处理 或 近似矩阵-向量积 ，而不是直接替代 \(A\)。例如，使用随机化方法构造一个预处理子 \(M \approx A^{-1}\)，然后求解 \(M A \mathbf{x} = M \mathbf{b}\)，此时矩阵 \(M A\) 的条件数改善，Chebyshev迭代收敛更快。 第五步：算法总结与复杂度分析 算法步骤 ： 随机化边界估计 ： 使用随机Lanczos（推荐）或随机幂迭代，运行 \(m\) 步，得到特征值边界估计 \(\tilde{a}, \tilde{b}\)。 加入安全边际：\(a = \tilde{a} - \epsilon, \quad b = \tilde{b} + \epsilon\)。 Chebyshev迭代求解 ： 初始化 \(\mathbf{x}_ 0, \mathbf{r}_ 0\)。 计算第一步：\(\mathbf{x}_ 1 = \mathbf{x}_ 0 + \frac{2}{b} \mathbf{r}_ 0\)。 对于 \(k = 1, 2, \dots\) 直到收敛： 计算 \(\omega_ {k+1}, \beta_ k\) 如上。 更新 \(\mathbf{x} {k+1}\) 和残差 \(\mathbf{r} {k+1}\)。 检查收敛。 计算复杂度 ： 边界估计：\(O(m \cdot \text{nnz}(A))\) 次浮点运算。 每轮迭代：一次矩阵-向量乘法（\(O(\text{nnz}(A))\)）和几次向量操作（\(O(n)\)）。 总成本：边界估计成本 + 迭代次数 × \(O(\text{nnz}(A))\)。 优势 ： 避免了精确特征值计算。 对于对称正定矩阵，收敛速度有保证（依赖于估计边界的质量）。 易于并行化（矩阵-向量乘法和向量操作均可并行）。 局限性 ： 主要适用于对称正定矩阵。对于非对称矩阵，需要转化为正规方程 \(A^\top A \mathbf{x} = A^\top \mathbf{b}\) 或使用其他多项式加速方法。 估计边界不准确时，可能收敛变慢甚至不稳定（安全边际可缓解）。 第六步：数值示例（概念性说明） 假设 \(A\) 是 \(1000 \times 1000\) 的对称正定稀疏矩阵，条件数 \(\kappa \approx 10^3\)。我们想解 \(A\mathbf{x} = \mathbf{b}\)。 运行随机Lanczos 30步，得到 \(\tilde{\lambda} {\min} \approx 0.95\)，\(\tilde{\lambda} {\max} \approx 950\)。 设置安全边际 \(\epsilon = 5\)，得到区间 \([ a, b] = [ 0.9, 955 ]\)。 运行Chebyshev迭代，参数 \(\omega_ k, \beta_ k\) 根据上述区间计算。 迭代约 \(k = 50\) 步后，残差下降至 \(10^{-6}\)（相对于初始残差）。 若使用经典CG方法，收敛步数约为 \(O(\sqrt{\kappa}) \approx 32\) 步。随机化Chebyshev可能略多几步（因边界估计误差），但避免了计算精确特征值，总体时间可能更优。 总结 随机化Chebyshev迭代法将 经典多项式加速 与 现代随机化技术 结合，为大型稀疏对称正定线性方程组提供了一种高效求解器。其核心贡献在于 通过随机化快速估计特征值边界 ，从而避免了传统Chebyshev迭代中对精确边界的需求。该方法在机器学习和科学计算中具有应用潜力，特别是当矩阵条件数较大且直接求解特征值不可行时。 值得注意的是，该算法可进一步与 随机采样 、 预处理技术 结合，以适应更广泛的矩阵类型和计算环境。