广义最小残差（GMRES）方法：原理与应用

玻尔百科

核心要点

GMRES通过迭代寻找一个近似解来求解线性系统，该近似解在一个被称为克雷洛夫子空间的扩展搜索空间上最小化残差范数。
实用的GMRES(m)算法每m步重启一次以控制内存成本，但这可能因丢弃有价值的信息而阻碍针对困难的非正规问题的收敛。
预处理和增广是通过将系统转换为一个更易于求解的系统，或在重启过程中保留关键信息，从而加速GMRES收敛的重要技术。
作为求解大型、稀疏、非对称系统的稳健求解器，GMRES在从计算流体动力学到控制理论和机器人学等领域中不可或缺。

引言

在科学与工程计算领域，许多复杂问题——从模拟流体流动到设计稳定的控制系统——最终都归结为求解形如 $Ax=b$ 的庞大线性方程组。当矩阵 $A$ 巨大、稀疏且非对称时，传统的直接方法便不再可行。这正是迭代方法的用武之地，其中功能最强大、用途最广泛的方法之一便是广义最小残差（GMRES）方法。GMRES提供了一种优雅而稳健的方式，通过在每一步做出最佳选择来找到一个日益精确的解。本文旨在不仅解释GMRES的工作原理，还要探讨它为何如此有效以及其局限性何在。

本次探索分为两大章节。首先，在“原理与机制”中，我们将深入GMRES的理论核心。我们将发现克雷洛夫子空间和最小多项式背后优雅的几何与代数思想，理解导致重启的GMRES( $m$ )变体的实际权衡，并学习它如何应对非正规矩阵的挑战性行为。随后，“应用与跨学科联系”章节将展示GMRES的实际应用。我们将了解它在计算流体动力学中的关键作用，其作为非线性求解器内部引擎的功能，以及它在机器人学和控制理论等领域出人意料的应用，同时强调预处理等技术不可或缺的作用。

原理与机制

既然我们已经了解了求解庞大线性方程组的巨大挑战，现在就让我们揭开广义最小残差方法（GMRES）的神秘面纱，欣赏其精巧的机制。你可能会想象，在一个拥有数百万维度的空间中找到一个精确解是一项无比复杂的任务。然而，GMRES背后的核心思想却异常简单，这证明了找到正确视角的力量。它始于一步，扩展为一个优美的几何结构，并最终达到一个强大的代数洞见。

克雷洛夫空间之旅

想象一下，你迷失在一片广阔的山区，目标是到达某个山谷的最低点。你有一个高度计，可以告诉你当前的“残差”——即你与目标点的垂直距离。一个明智的第一步是朝着最陡峭的下降方向移动。在我们的线性系统 $Ax=b$ 中，初始猜测 $x_0$ 的残差为 $r_0 = b - Ax_0$ 。这个残差向量指向一个误差方向。那么，为什么不朝那个方向迈出一步呢？

这正是GMRES的起点。它寻求一个更好的解 $x_1$ ，形式为 $x_1 = x_0 + \alpha r_0$ 。但是我们应该走多远？标量 $\alpha$ 的最佳值是什么？GMRES给出了一个清晰明确的答案：选择能使新残差 $r_1 = b - Ax_1$ 尽可能小的 $\alpha$ 。这个新残差的长度，即范数，是 $\|r_1\|_2 = \|r_0 - \alpha Ar_0\|_2$ 。我们只需找到最小化这个量的 $\alpha$ ——这是一个简单的微积分问题。

这个简单的想法已经揭示了一个深刻的真理。要让GMRES在第一步就找到精确解，残差必须变为零。这只有在 $r_0$ 和 $Ar_0$ 指向相同或相反方向时才可能，也就是说，当它们线性相关时。如果我们能找到一个 $\alpha$ 使得 $r_0 - \alpha Ar_0 = 0$ ，这意味着 $Ar_0 = (1/\alpha) r_0$ 。这正是特征向量的定义！因此，如果我们的初始残差 $r_0$ 恰好是矩阵 $A$ 的一个具有非零特征值 $\lambda$ 的特征向量，GMRES会巧妙地选择 $\alpha=1/\lambda$ 并在神奇的一步中收敛到精确解。

但如果 $r_0$ 不是特征向量呢？那么 $r_0$ 和 $Ar_0$ 指向不同的方向。我们能做的最好的事就是找到一个使 $r_1$ 尽可能短的 $\alpha$ ，但它不会是零。然而，请注意我们新残差的一个迷人之处： $r_1$ 是 $r_0$ 和 $Ar_0$ 的线性组合。那么，对于下一步，我们为什么要把自己局限于向新方向 $r_1$ 移动呢？我们已经有了两个有希望的方向， $r_0$ 和 $Ar_0$ 。为什么不在它们所张成的整个平面或子空间内搜索我们的下一次更新呢？

这便是GMRES的核心思想。在每一步 $k$ ，它不仅仅使用最新的残差，而是建立一个不断扩展的搜索方向库：

\mathcal{K}_k(A, r_0) = \mathrm{span}\{r_0, Ar_0, A^2r_0, \dots, A^{k-1}r_0\}

这个不断扩展的子空间被称为克雷洛夫子空间。可以把它想象成一个“回声”的画廊。初始残差 $r_0$ 是我们发出的声音。矩阵 $A$ 代表我们所处峡谷的复杂几何形状。 $Ar_0$ 是第一个回声， $A^2r_0$ 是回声的回声，依此类推。每一个新的回声都揭示了更多关于矩阵 $A$ 的结构。然后，GMRES通过考虑它迄今为止听到的所有回声中的所有信息，做出最优的移动。

在第 $k$ 步，GMRES在仿射空间 $x_0 + \mathcal{K}_k(A, r_0)$ 中找到“最佳”解 $x_k$ 。“最佳”始终意味着能最小化新残差范数 $\|b - Ax_k\|_2$ 的解。随着克雷洛夫子空间的增长，它在总解空间中占据了越来越大的部分，使GMRES有越来越大的机会找到一个残差非常小的点。

为了直观地看到这一点，考虑一个简单的 $2 \times 2$ 系统。非重启的GMRES保证最多在2步内找到精确解。第一步，它沿着由 $r_0$ 定义的直线搜索解。第二步，它在由 $\{r_0, Ar_0\}$ 张成的平面内搜索。由于整个空间只有二维，这个平面就是整个空间（假设 $r_0$ 和 $Ar_0$ 是独立的）。有了对整个空间的访问权限，它必然能找到精确解。对于一个 $N \times N$ 的系统，在最一般的情况下，克雷洛夫子空间 $\mathcal{K}_N(A, r_0)$ 将张成整个 $N$ 维空间，保证非重启的GMRES最多在 $N$ 次迭代内找到精确解。

最小多项式的魔力

在不断扩展的子空间中最小化向量长度的几何图像很美，但将视角转向代数，则揭示了该方法的真正天才之处。让我们看看残差 $r_k$ 的结构。由于我们的解是 $x_k = x_0 + z_k$ ，其中 $z_k \in \mathcal{K}_k(A, r_0)$ ，向量 $z_k$ 是克雷洛夫向量的线性组合： $z_k = c_0 r_0 + c_1 Ar_0 + \dots + c_{k-1}A^{k-1}r_0$ 。我们可以用矩阵 $A$ 的多项式更紧凑地表示这一点。如果我们定义一个多项式 $q_{k-1}(t) = c_0 + c_1 t + \dots + c_{k-1}t^{k-1}$ ，那么 $z_k = q_{k-1}(A)r_0$ 。

现在，让我们写出残差：

r_k = b - A x_k = b - A(x_0 + z_k) = (b-Ax_0) - Az_k = r_0 - A q_{k-1}(A) r_0

提出 $r_0$ ，我们得到：

r_k = (I - A q_{k-1}(A)) r_0

让我们定义一个新的多项式， $p_k(t) = 1 - t \cdot q_{k-1}(t)$ 。注意关于 $p_k(t)$ 的两件事：它的次数最多为 $k$ ，并且如果我们代入 $t=0$ ，我们得到 $p_k(0)=1$ 。我们的残差现在就简化为 $r_k = p_k(A)r_0$ 。

这便是GMRES的代数核心。在克雷洛夫子空间中找到最小化残差范数的向量的几何问题，与找到一个次数最多为 $k$ 且满足 $p_k(0)=1$ 的多项式 $p_k$ ，以最小化 $\|p_k(A)r_0\|_2$ 的范数是完全等价的。

这种多项式的观点非常强大。它告诉我们，GMRES不仅仅是在采取几何步骤；它是在隐式地构建一个最优的多项式算子来消灭初始误差。当我们可以找到这样一个多项式，使残差恰好为零时，就实现了向精确解的收敛。这在第 $k$ 步发生，如果我们能找到一个 $k$ 次多项式 $p_k$ 使得 $p_k(A)r_0 = 0$ 。能做到这一点的最低次多项式被称为 $A$ 关于 $r_0$ 的最小多项式。这个多项式的次数精确地告诉我们GMRES需要多少次迭代才能找到精确解。

这也解释了为什么非重启的GMRES必须在有限步内收敛。整个矩阵的最小多项式 $m_A(t)$ 具有 $m_A(A)$ 是零矩阵的性质。因此，对于任何 $r_0$ ， $m_A(A)r_0$ 始终为零。这个多项式的次数，比如说 $m$ ，是任何特定向量 $r_0$ 的最小多项式次数的上限。因此，GMRES总能在最多 $m$ 步内找到精确解，而 $m$ 总是小于或等于矩阵的大小 $N$ 。

内存成本与重启的艺术

所以，我们有了一个保证最多在 $N$ 步内奏效的方法。对于一个有一百万行的矩阵来说，这听起来像是一场灾难，而不是一个解决方案！但在实践中，我们希望GMRES远在 $N$ 步之前就给我们一个“足够好”的答案。更重要的是，我们所描述的“完全”GMRES存在一个令人望而生畏的实际问题。为了确保新的搜索方向是真正新的，GMRES使用一个称为Arnoldi过程的程序来为克雷洛夫子空间构建一个标准正交基。这意味着在第 $k$ 步，它必须存储 $k$ 个非常长的向量，每个向量的大小都是 $N$ 。

对于大规模问题，这是个致命弱点。如果 $N=10^6$ ，我们只运行 $k=300$ 次迭代，仅存储这些基向量就可能需要数GB的内存，远超许多计算机的内存容量。计算工作量也随每一步增长，因为每个新向量都必须与之前所有的向量正交。

实际的解决方案称为重启的GMRES，或GMRES( $m$ )。想法很简单：让GMRES运行一个固定的、适度的步数，比如 $m=50$ 。然后，取你找到的解 $x_m$ ，并以 $x_m$ 作为新的初始猜测，重新开始整个过程。这样可以保持内存和计算成本有界，因为你永远不需要存储超过 $m+1$ 个基向量。

但这种妥协是有代价的。通过重启，我们扔掉了克雷洛夫子空间及其包含的所有关于矩阵 $A$ 的宝贵信息。本质上，我们每 $m$ 步就让算法“失忆”一次。对于许多问题来说，这没问题。但对于一类特别棘手的矩阵，这可能是灾难性的。

驯服非正规巨兽：预处理与增广

GMRES的性能，特别是重启版本的性能，与矩阵 $A$ 的性质密切相关。如果一个矩阵是正规的（这一类别包括对称矩阵），它的特征值就足以说明一切。但许多来自现实世界物理问题的矩阵，例如来自流体动力学（对流-扩散问题）的矩阵，是非正规的。对于这些矩阵，特征值只是故事的一部分。

非正规矩阵可能表现出奇怪的瞬态行为，即在特征值所决定的长期行为生效之前，将矩阵应用于向量可能会暂时使其变得更大。这种行为与矩阵的特征向量非正交且可能近乎平行有关。衡量这一点的一个有用指标是特征向量矩阵的条件数 $\kappa(V)$ 。GMRES在此类矩阵上的收敛界限中包含这个 $\kappa(V)$ 项作为乘数。如果 $\kappa(V)$ 很大，即使特征值看起来漂亮地聚集在远离原点的地方，收敛也可能非常糟糕。

对于这些非正规问题，使用小 $m$ 重启GMRES尤其糟糕。GMRES( $m$ ) 被允许构建的短期多项式通常不够强大，无法越过初始的瞬态增长并开始减少误差。算法可能会停滞，从一次重启到下一次几乎没有进展。

能做些什么呢？这正是现代迭代方法艺术的用武之地。两种强大的策略是：

预处理： 其思想是“预处理”我们的系统，使其更容易求解。例如，使用右预处理，我们求解一个修改后的系统 $(AM^{-1})z = b$ ，然后通过 $x=M^{-1}z$ 恢复我们的解。矩阵 $M$ 是预处理器，被设计成 $A$ 的一个易于求逆的近似。如果 $M$ 是一个好的近似，新的矩阵 $AM^{-1}$ 会好得多——更接近单位矩阵——GMRES就能快速收敛。这种方法的一个巨大好处是，GMRES在每一步都最小化了我们原始问题的真实物理残差，因此我们的停止准则仍然有意义。
增广： 与其在重启时丢弃所有东西，不如我们识别出那些“制造麻烦”的方向——即GMRES难以处理的误差部分——并在重启时明确地将它们保留在我们的搜索空间中？这就是增广或降维GMRES背后的思想。通过在一个周期结束时分析克雷洛夫子空间，我们可以找到对应于接近零的、导致停滞的特征值的近似特征向量（通常称为谐波Ritz向量）。通过用这些向量增广下一个周期的搜索空间，我们为算法提供了克服非正规行为并高效解决这些困难问题所需的“长期记忆”。

归根结底，GMRES的故事是从一个简单、优雅的思想发展到一个复杂、实用的工具的历程。它结合了寻找最近点的几何直觉、多项式的代数力量，以及为使其在世界上最大的计算挑战中奏效所需的务实妥协。

应用与跨学科联系

我们花了一些时间来理解广义最小残差方法精巧的机制——它在不断扩展的搜索空间（克雷洛夫子空间）中寻找最佳解的优雅原理。但是，一台精美的机器只有在实际应用中才能真正被欣赏。这个算法到底能做什么？事实证明，答案是：几乎无所不能。

GMRES不仅仅是解决某个小众数学问题的专用工具。它是一把万能钥匙，开启了科学和工程领域中各种各样的问题。每当一个复杂系统——无论是流动的流体、振动的桥梁、金融市场，还是机器人的肢体——由方程描述时，通往其模拟的道路往往会导向一个庞大的线性方程组 $A x = b$ 。而矩阵 $A$ 常常是一个棘手、行为不佳的巨兽：巨大、稀疏、非对称，对更简单的方法来说简直是噩梦。这正是GMRES大放异彩的地方。让我们踏上一段旅程，看看GMRES帮助我们理解和控制的那些世界。

模拟的引擎：驯服自然法则

许多物理学的基本定律都以偏微分方程（PDEs）的形式表达，它们描述了热量、速度和压力等量如何随空间和时间变化。为了在计算机上求解这些方程，我们将空间和时间分割成精细的网格，并将PDE重写为一组巨大的耦合代数方程——我们熟悉的 $A x = b$ 。

不对称的挑战：对流与扩散

想象一下空气中的一缕烟。它既被风带着走（一个称为*对流的过程），又向四面八方扩散开来（一个称为扩散的过程）。这两种效应之间的平衡由一个称为佩克莱数（Péclet number）的无量纲量来捕捉。当扩散占主导地位时（低佩克莱数），问题是温和的，得到的矩阵 $A$ 通常是对称且行为良好的。但当对流占主导地位时——想象一下强劲而稳定的风——问题性质会发生巨大变化。矩阵变得高度非对称，或称非正规*。为对称系统设计的方法，如著名的共轭梯度法，会彻底失败。

然而，GMRES正是为此而生。它的收敛不依赖于对称性。通过直接最小化残差，即使面对由对流主导的物理现象所引起的强非正规性，它也能有条不紊地逐步求解。研究当我们提高佩克莱数时GMRES的收敛性如何变化，为我们提供了一个直接窗口，来观察该方法的威力及其局限性，这是对任何旨在模拟现实世界的工具的关键压力测试。

流体的非定常世界

考虑模拟水流过管道或空气流过机翼的任务。对于不可压缩流体，我们必须求解斯托克斯（Stokes）或纳维-斯托克斯（Navier-Stokes）方程。当离散化后，这些方程通常会导致所谓的*鞍点系统*。这些问题的矩阵通常是不定的，意味着它们既有正特征值也有负特征值。这是另一个许多迭代方法会迷失的领域。GMRES对特征值的符号漠不关心，只关心最小化残差，因此能够驾驭这片险恶的地形。它已成为计算流体动力学（CFD）中不可或缺的工具，帮助我们设计更高效的飞机并理解复杂的天气模式。

不可或缺的伙伴：预处理

即使是像GMRES这样强大的引擎，如果问题过于困难，也可能力不从心。一个来自真实世界模拟的原始、未预处理的系统可能需要如此多的迭代才能求解，以至于在计算上变得不可行。解决方案不是放弃GMRES，而是给它一个伙伴：预处理器。

预处理器 $M$ 是我们困难矩阵 $A$ 的一个近似，关键在于它易于求逆。我们可能不是求解 $A x = b$ ，而是求解等价的系统 $A M^{-1} y = b$ ，然后通过 $x = M^{-1} y$ 恢复我们的解。这被称为*右预处理*。目标是使新矩阵 $A M^{-1}$ 比原始的 $A$ “更好”——其特征值更聚集，非正规性更低——从而使GMRES能以少得多的步数收敛。

一个主力预处理策略是不完全LU（ILU）分解。对稀疏矩阵进行完全LU分解可能会产生灾难性的密集结果且成本高昂。ILU分解通过在计算过程中有策略地丢弃一些元素来计算近似的L和U，从而保持稀疏性。即使是这种“不完全”的近似，也常常足以将GMRES的速度提高几个数量级。
在一个算法巧思的美妙转折中，我们甚至可以用一种迭代方法来预处理另一种。几次简单方法（如逐次超松弛法，SOR）的扫描可以起到应用近似逆的作用，从而成为外部GMRES循环的有效预处理器。这显示了这些数值工具非凡的模块化和灵活性。

预处理的选择是一门艺术，但它由深刻的原则指导。例如，右预处理有一个很好的特性，即GMRES最小化的残差是原始问题的真实残差。左预处理（ $M^{-1} A x = M^{-1} b$ ）有时更容易分析，但它最小化的是一个“预处理过的”残差，要得到真实残差需要额外的一步。理解这些细微差别是实际、有效使用GMRES的关键。

超越网格：通用的问题求解器

GMRES的影响力远远超出了PDE模拟的结构化网格。对于任何可以被线性化的问题，它都是一个基础工具。

非线性求解器的核心

大多数现实世界的问题都是非线性的。为了求解非线性系统 $F(x) = 0$ ，我们拥有的最强大的方法是牛顿法。这是一个迭代过程：从当前猜测 $x_k$ 出发，我们用切线（或切平面）来近似非线性函数，并找到该切线与零点的交点。这给了我们下一个更好的猜测 $x_{k+1}$ 。这一步的计算需要求解一个线性系统： $J(x_k) s_k = -F(x_k)$ ，其中 $J(x_k)$ 是雅可比矩阵。

对于大问题，在每个牛顿步骤中都精确求解这个线性系统是极其昂贵的。这就是*非精确牛顿法的用武之地。我们可以使用GMRES来近似地*求解线性系统——只需足够精确，以便在外部的非线性问题上取得良好进展。通过仔细控制内部GMRES求解的容差，我们可以实现惊人的效率，将一个棘手的非线性问题转化为一系列可控的线性问题。这使得GMRES成为科学计算和优化庞大机器中的一个关键引擎，使得解决以前无法企及的问题成为可能。

控制复杂系统与“无矩阵”的力量

工程师如何确保一个复杂系统——一架飞机、一个电网、一个化学反应器——是稳定的？控制理论的基石之一是李雅普诺夫方程，一个形如 $A X + X A^{T} = -C$ 的矩阵方程。求解这个方程以得到矩阵 $X$ ，可以告诉我们由矩阵 $A$ 控制的系统的稳定性。

对于一个有 $n$ 个状态的系统，这个矩阵方程可以转化为一个大小为 $n^2 \times n^2$ 的线性系统。如果 $n=1000$ ，这个矩阵将有一万亿个元素！存储这样一个矩阵是不可能的。在这里，我们见证了GMRES最深刻的特性之一：它可以以无矩阵（matrix-free）的方式实现。GMRES不需要看到矩阵；它只需要一个“黑匣子”函数，给定一个向量 $v$ ，返回乘积 $Av$ 。对于李雅普诺夫方程，我们可以编写这样一个函数，而无需构建那个巨大的 $n^2 \times n^2$ 算子。GMRES通过与这个黑匣子交互，可以像矩阵明确存在一样求解该系统。这种无矩阵的理念是一次范式转变，使我们能够处理规模巨大的问题。

机器人学与运动几何学

GMRES的多功能性甚至延伸到了机器人学、计算机视觉和图形学领域。想象一个机械臂试图将一个工具与目标对齐。这可以被构建为寻找一个小的校正旋转。使用优雅的四元数代数来表示旋转，这个校正问题可以被线性化为一个小的方程组， $A \delta\omega = b$ ，其中 $\delta\omega$ 是我们想要找到的小旋转向量。这个系统可能是病态的，但通过用GMRES求解一个稍作修改的、正则化的版本，我们可以找到一个稳定而准确的解。这是一个美丽的示范，展示了抽象的代数结构（四元数）如何导致一个非常适合我们最小残差英雄的 конкрет 数值问题。

求解器的艺术：实践智慧与未来展望

最后，将GMRES置于其适当的背景中非常重要。它是一个强大的工具，但不是唯一的工具。而且，像任何工具一样，它也在不断被改进和调整以适应新的挑战。

求解器的世界：GMRES与BiCGSTAB的对决

对于非对称系统，另一个流行的方法是双共轭梯度稳定法（BiCGSTAB）。这两种方法有不同的理念。GMRES是谨慎的优化者：在一个周期内的每一步，它都保证在其搜索空间内得到最佳解（在最小残差的意义上）。这种最优性是有代价的：计算工作量和内存随着周期中的每一步而增长。为了控制这一点，我们每 $m$ 步重启GMRES（使用GMRES( $m$ )），但这意味着我们丢弃了有价值的信息，如果 $m$ 太小，这可能会减慢甚至导致收敛停滞。

另一方面，BiCGSTAB是一个灵活的机会主义者。它使用短的、固定成本的递推，并且不要求每一步都达到最优。它的收敛过程可能看起来不规则，残差范数会上下跳动。然而，因为它不重启，它会不断地利用其过去的所有工作。在好的预处理器使问题变得“容易”的情况下，或者在GMRES( $m$ )因重启值太小而停滞的情况下，BiCGSTAB低廉且恒定的每次迭代成本可能使其更快地冲向终点。选择正确的求解器是一门实践艺术，需要根据手头问题的具体结构来决定。

超级计算机时代的GMRES

随着我们将计算的前沿推向大规模并行超级计算机，一个新的挑战出现了。处理器之间通信所需的时间（延迟）可能会远远超过它们进行实际计算的时间。一个需要频繁进行全局通信的算法将会很慢，无论处理器有多快。

经典的GMRES，由于其每一步都有内积计算，是通信密集型的。这催生了避免通信算法（如CA-GMRES）的发展。其思想是重构算法，在每个处理器上本地执行 $s$ 步的计算，然后执行一次单一的、更大规模的通信。这用额外的浮点运算换取了延迟成本的大幅降低。通过分析这种权衡，甚至可以推导出在给定并行架构上最小化总时间的最优块大小 $s$ 。这项工作表明，GMRES的故事尚未结束。该算法在不断演进，将其优雅的核心原理适应于高性能计算不断变化的格局。

从其在克雷洛夫子空间中最小化残差的数学核心出发，我们已经看到GMRES的影响力扩展到计算世界的几乎每一个角落。它的美不仅在于其优雅的理论，更在于其作为一种稳健、通用且必不可少的科学发现工具，所具有的深远且不断扩展的实用性。