GMRES

Intuition Beginner

Suppose you face a giant system of equations $Ax=b$, so large you can never write the matrix down, but you can do one cheap thing: hand the matrix a vector and get back the matrix times that vector. From this one power you grow a small space of candidate answers — the Krylov space of 43.07.02 — built from your starting guess and its successive products by the matrix. The question of this unit is: among all the candidate answers in that growing space, which one is best?

GMRES answers with a precise rule. For any candidate answer $x$, the residual $b-Ax$ measures how badly the equation is missed: a residual of zero means $x$ is exact. GMRES picks, at each stage, the candidate whose residual is as short as possible. That is why it is called the generalised minimal residual method. As the space grows by one direction each stage, the best residual can only get shorter or stay the same — never longer. So the method marches steadily toward the answer.

There is a clever bookkeeping trick that makes this cheap. Building the Krylov space the tidy way (the Arnoldi process) turns the giant minimisation into a tiny one: at stage $m$ you only solve a small lopsided system with $m+1$ rows and $m$ columns. And you can read off how short the current residual is without ever assembling the full candidate answer — a single number falls out of the bookkeeping for free, so you know when to stop before paying to form the solution.

GMRES works for a general matrix, even a lopsided non-symmetric one. When the matrix is symmetric and positive-definite, a cheaper cousin called conjugate gradients does the same job with short-term memory; GMRES is the price you pay for generality.

Visual Beginner

The picture is a growing space of candidate answers, with GMRES standing at the true target and, at each stage, picking the candidate closest to it — so the leftover miss shrinks step by step.

Read the table top to bottom. Each row is one stage. The middle column is the length of the shortest residual achievable from the space built so far; it never grows. The right column notes the size of the tiny least-squares problem solved at that stage — always one row taller than it is wide.

stage $m$	shortest residual length	small problem solved
1	$|r_1|\le |r_0|$	$2\times 1$
2	$|r_2|\le |r_1|$	$3\times 2$
3	$|r_3|\le |r_2|$	$4\times 3$
4	$|r_4|\le |r_3|$	$5\times 4$

The takeaway: GMRES is "find the nearest point of a growing space to a fixed target". Because each new space contains the last, the nearest point can only get nearer. The work at each stage is a small, tall least-squares fit, and the size of the leftover at the bottom of that fit is exactly the residual length you are driving to zero.

Worked example Beginner

Take the small system $$ A = \begin{pmatrix} 2 & 1 \ 0 & 2 \end{pmatrix}, \qquad b = \begin{pmatrix} 1 \ 1 \end{pmatrix}, $$ and start from the guess $x_0=(0,0)$, so the starting residual is $r_0=b-A{x}_0=b=(1,1)$, of length $\sqrt{2}$.

Stage 1 (the first Krylov direction). The first tidied direction is $r_0$ rescaled to length one: $q_1=(1,1)/\sqrt{2}$. GMRES looks for the best candidate of the form $x_1=x_0+\tau q_1$ — a single number $\tau$ to choose. We want to make $b-A{x}_1$ as short as possible.

Compute $A{q}_1=\frac{1}{\sqrt{2}}(2+1,0+2)=\frac{1}{\sqrt{2}}(3,2)$. The residual of the candidate is $r_0-\tau A{q}_1=(1,1)-\frac{\tau }{\sqrt{2}}(3,2)$.

Choosing $\tau$ to minimise the length of $(1-3\tau /\sqrt{2},1-2\tau /\sqrt{2})$ is a one-variable least-squares fit. Setting the derivative of the squared length to zero: the squared length is $(1-3s{)}^2+(1-2s{)}^2$ with $s=\tau /\sqrt{2}$, whose derivative is $-6(1-3s)-4(1-2s)=-10+26s$, zero at $s=10/26=5/13$. The shortest residual then has squared length $(1-15/13{)}^2+(1-10/13{)}^2=(2/13{)}^2+(3/13{)}^2=13/169=1/13$, so $\Vert r_1\Vert =1/\sqrt{13}\approx 0.277$.

What this tells us: one stage cut the residual from $\sqrt{2}\approx 1.414$ down to about $0.277$, just by choosing the single best step along the first Krylov direction. A second stage would add the next direction and drive the residual to zero, since the space then fills all of the plane. Each stage solved only a tiny minimisation, never touching $A$ except through one matrix-times-vector product.

Check your understanding Beginner

Formal definition Intermediate+

Let $A\in {\mathbb{F}}^{n\times n}$ with $\mathbb{F}\in \{\mathbb{R},\mathbb{C}\}$ be nonsingular and possibly nonsymmetric, let $b\in {\mathbb{F}}^n$, and fix an initial guess $x_0$ with residual $r_0=b-A{x}_0\ne 0$. The $m$-th GMRES iterate is the unique minimiser $$ x_m = \operatorname*{arg,min}_{x \in x_0 + \mathcal{K}_m(A, r_0)} |b - A x|_2, \qquad \mathcal{K}_m(A, r_0) = \operatorname{span}{r_0, A r_0, \dots, A^{m-1} r_0}, $$ where ${\mathcal{K}}_m(A,r_0)$ is the Krylov subspace of 43.07.02 generated by the residual. Equivalently, writing $x=x_0+z$ with $z\in {\mathcal{K}}_m$, GMRES minimises $\Vert r_0-Az{\Vert }_2$ over $z\in {\mathcal{K}}_m$ ^{[Saad, Y. & Schultz, M. H. — GMRES: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems]}.

The Arnoldi reduction. Run the Arnoldi process of 43.07.02 on $A$ with starting vector $r_0$. It produces, before breakdown, an orthonormal basis $Q_m=[q_1\cdots q_m]$ of ${\mathcal{K}}_m$ with $q_1=r_0/\beta$, $\beta =\Vert r_0{\Vert }_2$, together with the $(m+1)\times m$ upper-Hessenberg matrix ${\tilde{H}}_m$ satisfying the Arnoldi relation $A{Q}_m=Q_{m+1}{\tilde{H}}_m$. Parametrise $z=Q_{m}y$ with $y\in {\mathbb{F}}^m$. Then $$ r_0 - A z = \beta q_1 - A Q_m y = \beta q_1 - Q_{m+1}\tilde H_m y = Q_{m+1}\big(\beta e_1 - \tilde H_m y\big), $$ using $q_1=Q_{m+1}{e}_1$. Since $Q_{m+1}$ has orthonormal columns it preserves the Euclidean norm, so $$ |b - A x|_2 = |r_0 - A z|_2 = \big|,\beta e_1 - \tilde H_m y,\big|2. $$ The large minimisation collapses to the small $(m+1)\times m$ least-squares problem $\min{y}|\beta e_1 - \tilde H_m y|_2$,ofthekindtreatedin[\mathrm{43.04.01}].Theminimiser$y_m$gives$x_m = x_0 + Q_m y_m$,andtheresidualnormequalstheoptimalleast-squaresresidual$|\beta e_1 - \tilde H_m y_m|_2$.

Givens-rotation solution and the cheap residual. Because ${\tilde{H}}_m$ is upper-Hessenberg, its QR factorisation is built by $m$ plane (Givens) rotations $G_1,\dots ,G_m$, each zeroing one subdiagonal entry $h_{j+1,j}$ — the incremental QR of 43.04.01. Applying the accumulated rotations to $\beta e_1$ yields a vector $g\in {\mathbb{F}}^{m+1}$ whose first $m$ entries feed a triangular back-substitution for $y_m$ and whose last entry $g_{m+1}$ satisfies $$ |b - A x_m|2 = |g{m+1}|. $$ The residual norm is therefore available at every stage without forming $x_m$: one monitors $|g_{m+1}|$ and only assembles $x_m=x_0+Q_m{y}_m$ once a stopping tolerance is met. Each new stage appends one rotation, so the update is incremental.

Restarted GMRES. Storing $Q_m$ costs $O(mn)$ and the orthogonalisation cost grows like $O(m^{2}n)$; for large $m$ this is prohibitive. Restarted GMRES, written GMRES$(k)$, runs $k$ stages, forms $x_k$, discards $Q_k$, and restarts the whole process with $x_0\leftarrow x_k$ and $r_0\leftarrow b-A{x}_k$. This caps storage at $O(kn)$ at the cost of losing the global optimality across restart boundaries.

Counterexamples to common slips

• GMRES minimises the residual, not the error. The quantity made small is $\Vert b-A{x}_m{\Vert }_2=\Vert A(x_{\star }-x_m){\Vert }_2$, not $\Vert x_{\star }-x_m{\Vert }_2$. For an ill-conditioned $A$ a small residual can coexist with a large error; the error bound carries a factor $\kappa (A)$.
• Breakdown is a feature. The Arnoldi subdiagonal vanishing, $h_{m+1,m}=0$, is a lucky breakdown: it signals that ${\mathcal{K}}_m$ is $A$-invariant and that $x_m$ is the exact solution, residual zero. Unlike the nonsymmetric Lanczos process, GMRES cannot suffer a serious breakdown.
• The spectrum alone need not predict convergence. For a non-normal $A$, eigenvalues clustered away from the origin do not guarantee fast GMRES convergence; the eigenvector conditioning $\kappa (V)$, the field of values, or the pseudospectrum enters. The diagonalisable bound is only a bound, and a loose one when $V$ is ill-conditioned.
• Restarted GMRES$(k)$ can stagnate. Plain (full) GMRES converges in at most $n$ steps in exact arithmetic, but GMRES$(k)$ can stall short of the solution, the residual flat across many restart cycles, when the discarded directions were the ones that mattered.

Key theorem with proof Intermediate+

The signature result is that the GMRES minimisation over the Krylov subspace reduces exactly to the small Hessenberg least-squares problem, that the residual norm is the last rotated component, and that the residual is monotonically non-increasing with finite termination.

Theorem (GMRES least-squares reduction, monotonicity, and finite termination). Let $A$ be nonsingular, $r_0=b-A{x}_0\ne 0$, $\beta =\Vert r_0{\Vert }_2$, and run Arnoldi on $(A,r_0)$ producing $Q_{m+1}$ and the unreduced $(m+1)\times m$ Hessenberg ${\tilde{H}}_m$. Then the GMRES iterate satisfies $$ x_m = x_0 + Q_m y_m, \qquad y_m = \operatorname*{arg,min}_{y \in \mathbb{F}^m}\big|\beta e_1 - \tilde H_m y\big|_2, \qquad |b - A x_m|2 = \min_y\big|\beta e_1 - \tilde H_m y\big|2. $$ *The residual norms are non-increasing, $|b - A x{m+1}|2 \le |b - A x_m|2$,andthereisaleastindex$\nu \le n$(thegradeof$r_0$)atwhich$h{\nu+1,\nu} = 0$and$x\nu = x\star$ exactly* ^{[Trefethen, L. N. & Bau, D. — Numerical Linear Algebra (SIAM, 1997)]}.

Proof. The reduction is the computation of the Formal definition. Every $x\in x_0+{\mathcal{K}}_m$ is $x=x_0+Q_{m}y$ for a unique $y\in {\mathbb{F}}^m$ since $Q_m$ has full column rank. The Arnoldi relation $A{Q}_m=Q_{m+1}{\tilde{H}}_m$ and $r_0=\beta q_1=\beta Q_{m+1}{e}_1$ give $$ b - A x = r_0 - A Q_m y = Q_{m+1}(\beta e_1 - \tilde H_m y). $$ The columns of $Q_{m+1}$ are orthonormal, so $\Vert Q_{m+1}u{\Vert }_2=\Vert u{\Vert }_2$ for every $u\in {\mathbb{F}}^{m+1}$; hence $\Vert b-Ax{\Vert }_2=\Vert \beta e_1-{\tilde{H}}_{m}y{\Vert }_2$. Minimising the left side over $x\in x_0+{\mathcal{K}}_m$ is therefore the same as minimising the right side over $y\in {\mathbb{F}}^m$, which is the stated least-squares problem; its minimiser $y_m$ yields $x_m=x_0+Q_m{y}_m$ and residual norm equal to the optimal value.

Monotonicity follows from the nesting ${\mathcal{K}}_m\subseteq {\mathcal{K}}_{m+1}$. The feasible set $x_0+{\mathcal{K}}_m$ for stage $m$ is contained in that for stage $m+1$, so the minimum of $\Vert b-Ax{\Vert }_2$ over the larger set is at most the minimum over the smaller: $\Vert b-A{x}_{m+1}{\Vert }_2\le \Vert b-A{x}_m{\Vert }_2$.

For finite termination, the Krylov dimensions $\dim {\mathcal{K}}_m$ strictly increase by one until the grade $\nu$ of $r_0$ — the least $m$ with ${\mathcal{K}}_m={\mathcal{K}}_{m+1}$, established in 43.07.02 — and at that step the Arnoldi remainder vanishes, $h_{\nu +1,\nu }=0$, so ${\tilde{H}}_{\nu }=H_{\nu }$ is the square $\nu \times \nu$ Hessenberg with $A{Q}_{\nu }=Q_{\nu }{H}_{\nu }$. The least-squares problem becomes the square system $H_{\nu }y=\beta e_1$. Since $A$ is nonsingular and ${\mathcal{K}}_{\nu }$ is $A$-invariant, the restriction of $A$ to ${\mathcal{K}}_{\nu }$ is nonsingular, so $H_{\nu }=Q_{\nu }^{\ast }A{Q}_{\nu }$ is nonsingular and the system has a solution $y_{\nu }$ with zero residual. Then $b-A{x}_{\nu }=Q_{\nu +1}(\beta e_1-{\tilde{H}}_{\nu }{y}_{\nu })=0$, so $x_{\nu }=x_{\star }$, and $\nu \le n$ because $\dim {\mathcal{K}}_{\nu }\le n$. $\square$

Bridge. This reduction is the foundational reason GMRES is practical at all: it builds toward the convergence theory of the Advanced results by replacing an $n$-dimensional minimisation with an $(m+1)\times m$ Hessenberg least-squares problem whose answer carries its own error estimate, the last rotated Givens component. This is exactly the Arnoldi factorisation of 43.07.02 read as a solver rather than an eigenvalue probe — there the projected $H_m$ supplies Ritz values, here the augmented ${\tilde{H}}_m$ supplies the residual-minimising step — and the monotone residual decrease generalises the monotone Rayleigh-quotient bound of that unit from extremal eigenvalues to the full residual. The least-squares inner solve is exactly the incremental Givens QR of 43.04.01 applied to a Hessenberg matrix, so the per-step cost is one rotation. The construction appears again in 43.07.04, where for symmetric positive-definite $A$ the same Krylov optimality is achieved with a short three-term recurrence instead of a growing orthonormal basis; putting these together, the central insight is that GMRES is "least squares on the Arnoldi Hessenberg", and the bridge is that every quantity a practitioner needs — the iterate, the residual norm, the stopping decision — is read off this one small projected problem one column at a time.

Exercises Intermediate+

Advanced results Master

Theorem 1 (minimal-residual optimality and the residual polynomial). Among all $x\in x_0+{\mathcal{K}}_m(A,r_0)$ the GMRES iterate $x_m$ is the unique minimiser of $\Vert b-Ax{\Vert }_2$, and its residual obeys the optimality characterisation $$ |r_m|2 = \min{p \in P_m,, p(0)=1}|p(A) r_0|_2, $$ where $P_m$ is the space of polynomials of degree at most $m$. Uniqueness is convexity: $y\mapsto \Vert \beta e_1-{\tilde{H}}_{m}y{\Vert }_2$ is a strictly convex function on the range of the full-column-rank ${\tilde{H}}_m$ (an unreduced Hessenberg has rank $m$), so its minimiser is unique, and the affine map $y\mapsto x_0+Q_{m}y$ transports uniqueness to $x_m$. The polynomial form follows because the GMRES residual is $p(A)r_0$ with $p(0)=1$, and minimising over the subspace is minimising over such polynomials. This recasts GMRES convergence as a problem in approximation theory on $A$: how well can a degree-$m$ polynomial, pinned to value one at the origin, annihilate $r_0$ ^{[Saad, Y. & Schultz, M. H. — GMRES: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems]}.

Theorem 2 (diagonalisable convergence bound and the clustering heuristic). If $A=V\Lambda V^{-1}$ is diagonalisable, then $$ \frac{|r_m|2}{|r_0|2} \le \kappa_2(V)\ \min{p\in P_m,, p(0)=1}\ \max{\lambda \in \sigma(A)}|p(\lambda)|. $$ The minimax factor is small whenever the spectrum splits into a few tight clusters bounded away from the origin: a polynomial with one root near each cluster of total degree equal to the number of clusters is already tiny on $\sigma (A)$, so GMRES converges in roughly as many stages as there are clusters. The factor ${\kappa }_2(V)$ quantifies the eigenvector conditioning; for normal $A$ it is one and the bound is essentially tight, but for non-normal $A$ it can be astronomically large, and then the bound says little. The clustering heuristic — "GMRES converges fast when the eigenvalues cluster" — is reliable for normal or mildly non-normal matrices and is exactly where preconditioning aims, to reshape $\sigma (A)$ into clusters ^{[Saad, Y. — Iterative Methods for Sparse Linear Systems (2nd ed.)]}.

Theorem 3 (field-of-values bound for the non-normal case). Let $W(A) = {x^ A x : |x|_2 = 1}$bethefieldofvalues(numericalrange)of$A$,andsuppose$0 \notin W(A)$,sothedistance$\nu = \operatorname{dist}(0, W(A)) > 0$and$W(A)$liesinadiscofradius$R$aboutsomecentre.ThenGMRESconvergeslinearlywitharatedeterminedby$W(A)$ratherthanby$\sigma(A)$ alone; a representative Elman-type bound is* $$ \frac{|r_m|_2}{|r_0|_2} \le \Big(1 - \frac{\nu^2}{|A|2^2}\Big)^{m/2}. $$ When $A$ is far from normal the spectrum can be deeply misleading — eigenvalues bunched away from zero while the field of values straddles the origin, in which case GMRES stagnates despite a "good" spectrum. The field of values, and more finely the pseudospectrum $\sigma\varepsilon(A) = {z : |(zI - A)^{-1}|_2 \ge \varepsilon^{-1}}$, are the correct convergence-controlling sets: a min-max polynomial small on a pseudospectral region enclosing the spectrum at a safe distance from the origin yields a genuine residual bound by a contour-integral (Cauchy) estimate ^{[Greenbaum, A. — Iterative Methods for Solving Linear Systems (SIAM, 1997)]}.

Theorem 4 (any-curve theorem; the spectrum does not determine the GMRES curve). Given any set of $n$ nonzero complex numbers (a prescribed spectrum) and any non-increasing sequence $\Vert r_0\Vert \ge \Vert r_1\Vert \ge \cdots \ge \Vert r_{n-1}\Vert >\Vert r_n\Vert =0$, there exists a matrix $A$ with that spectrum and a right-hand side $b$ such that GMRES applied to $Ax=b$ produces exactly that residual sequence. The construction (Greenbaum, Pták, Strakoš) decouples the eigenvalues from the convergence history entirely: any admissible residual curve is compatible with any prescribed spectrum, so no theorem bounding GMRES convergence by the eigenvalues alone can exist for general matrices. This is the sharp statement of why the diagonalisable bound's ${\kappa }_2(V)$ prefactor cannot be removed and why non-normal convergence analysis must invoke the field of values or pseudospectra. For the symmetric case the situation reverts: for $A=A^{\ast }$ the eigenvector matrix is unitary, the prefactor is one, and the spectrum does determine the convergence curve ^{[Greenbaum, A. — Iterative Methods for Solving Linear Systems (SIAM, 1997)]}.

Synthesis. GMRES is one construction — minimise $\Vert b-Ax{\Vert }_2$ over the affine Krylov space $x_0+{\mathcal{K}}_m(A,r_0)$ — viewed under one invariant, the residual polynomial $p(A)r_0$ with $p(0)=1$, and the foundational reason the method works is that this minimisation is, by the Arnoldi factorisation of 43.07.02, a small Hessenberg least-squares problem whose residual norm is read off a single Givens component. This is exactly the Arnoldi eigenvalue engine of that unit turned to a different end: there the projected $H_m$ delivers Ritz values approximating the spectrum, here the augmented ${\tilde{H}}_m$ delivers the residual-minimising step, and the same nesting ${\mathcal{K}}_m\subseteq {\mathcal{K}}_{m+1}$ that made the extremal Ritz values monotone makes the GMRES residual monotone. The polynomial recasting of Theorem 1 generalises the static eigenvalue-approximation picture into a dynamical approximation-theory one: the convergence rate is a min-max polynomial value on a spectral set, and the central insight is that whether that set is the spectrum, the field of values, or the pseudospectrum is decided by the normality of $A$.

The non-normal subtlety is where the bridge to honest practice is built: Theorem 2's clustering heuristic is reliable only through the prefactor ${\kappa }_2(V)$, Theorem 3 replaces the spectrum by the field of values when $A$ is far from normal, and Theorem 4 shows the eigenvalues alone determine nothing about the GMRES curve for general matrices — the any-curve construction is dual to the diagonalisable bound's looseness, each expressing that $p(A)$ is not $p(\Lambda )$. Putting these together, GMRES is the general-matrix residual-minimiser whose symmetric-positive-definite specialisation reappears in 43.07.04 as conjugate gradients with a short recurrence and a spectrum-determined error bound; the bridge is that the stationary-iteration splitting $M$ of 43.07.01, powerless on its own against a spectral radius creeping to one, returns here as a preconditioner whose only job is to cluster the spectrum or shrink the field of values so the GMRES residual polynomial can do its work in a handful of stages.

Full proof set Master

Proposition 1 (least-squares reduction and norm identity). With $Q_{m+1},{\tilde{H}}_m$ from Arnoldi on $(A,r_0)$ and $\beta =\Vert r_0{\Vert }_2$, the GMRES iterate is $x_m=x_0+Q_m{y}_m$ with $y_m=\arg {\min }_y\Vert \beta e_1-{\tilde{H}}_{m}y{\Vert }_2$, and $\Vert b-A{x}_m{\Vert }_2=\Vert \beta e_1-{\tilde{H}}_m{y}_m{\Vert }_2$.

Proof. Each $x\in x_0+{\mathcal{K}}_m$ is $x_0+Q_{m}y$ for a unique $y$ because $Q_m$ has orthonormal (hence independent) columns spanning ${\mathcal{K}}_m$. The Arnoldi relation $A{Q}_m=Q_{m+1}{\tilde{H}}_m$ and $r_0=\beta q_1=\beta Q_{m+1}{e}_1$ give $b-Ax=r_0-A{Q}_{m}y=Q_{m+1}(\beta e_1-{\tilde{H}}_{m}y)$. Orthonormality of $Q_{m+1}$ yields $\Vert b-Ax{\Vert }_2=\Vert \beta e_1-{\tilde{H}}_{m}y{\Vert }_2$. Minimising both sides over their (corresponding) feasible sets is the same problem; the minimiser $y_m$ gives $x_m$ and the residual identity. $\square$

Proposition 2 (monotonicity of the residual). $\Vert r_{m+1}{\Vert }_2\le \Vert r_m{\Vert }_2$ for all $m$ before termination.

Proof. The feasible affine spaces are nested, $x_0+{\mathcal{K}}_m\subseteq x_0+{\mathcal{K}}_{m+1}$, since ${\mathcal{K}}_m\subseteq {\mathcal{K}}_{m+1}$. The minimum of a fixed function over a larger set cannot exceed its minimum over a subset, so $\Vert r_{m+1}{\Vert }_2={\min }_{x\in x_0+{\mathcal{K}}_{m+1}}\Vert b-Ax{\Vert }_2\le {\min }_{x\in x_0+{\mathcal{K}}_m}\Vert b-Ax{\Vert }_2=\Vert r_m{\Vert }_2$. $\square$

Proposition 3 (lucky breakdown equals exact solution). If $h_{m+1,m}=0$ at stage $m$ and $A$ is nonsingular, then $x_m=x_{\star }$ with $r_m=0$, and this is the only way GMRES can break down.

Proof. Breakdown closes the Arnoldi relation to $A{Q}_m=Q_m{H}_m$, so ${\mathcal{K}}_m=\mathrm{range}(Q_m)$ is $A$-invariant. Restricted to the $m$-dimensional invariant ${\mathcal{K}}_m$, the nonsingular $A$ acts bijectively, so $H_m=Q_m^{\ast }A{Q}_m$ is nonsingular and $H_m{y}_m=\beta e_1$ has the exact solution $y_m=\beta H_m^{-1}{e}_1$. The residual is $b-A{x}_m=Q_m(\beta e_1-H_m{y}_m)=0$. No serious breakdown is possible: the Arnoldi step fails only by producing the zero remainder $w=0$, which is precisely $h_{m+1,m}=0$, and that case has just been shown to be convergence. (Contrast the two-sided nonsymmetric Lanczos process, where a serious breakdown — a vanishing bilinear product with nonzero vectors — can halt progress without convergence.) $\square$

Proposition 4 (finite termination at the grade). Let $\nu$ be the grade of $r_0$ with respect to $A$. Then GMRES produces the exact solution at stage $\nu$, and $\nu \le \mathrm{deg}{m}_A\le n$, with $\mathrm{deg}{m}_A$ the degree of the minimal polynomial of $A$.

Proof. By 43.07.02, the Krylov dimensions strictly increase until the grade $\nu$, where ${\mathcal{K}}_{\nu }={\mathcal{K}}_{\nu +1}$ is $A$-invariant and $h_{\nu +1,\nu }=0$. Proposition 3 then gives $x_{\nu }=x_{\star }$. The grade is the degree of the minimal polynomial of $A$ restricted to the cyclic subspace generated by $r_0$, which divides $m_A$, so $\nu \le \mathrm{deg}{m}_A\le n$. In exact arithmetic GMRES is therefore a finite method; its value is that the residual usually falls below tolerance long before stage $\nu$. $\square$

Proposition 5 (residual-polynomial optimality and the diagonalisable bound). The GMRES residual satisfies $\Vert r_m{\Vert }_2={\min }_{p\in P_m,p(0)=1}\Vert p(A)r_0{\Vert }_2$, and for diagonalisable $A=V\Lambda V^{-1}$, $\frac{\Vert r_m{\Vert }_2}{\Vert r_0{\Vert }_2}\le {\kappa }_2(V){\min }_{p\in P_m,p(0)=1}{\max }_{\lambda \in \sigma (A)}|p(\lambda )|$.

Proof. The Krylov subspace is ${\mathcal{K}}_m=\{q(A)r_0:\mathrm{deg}q<m\}$, so any $x\in x_0+{\mathcal{K}}_m$ has residual $b-Ax=(I-Aq(A))r_0=p(A)r_0$ with $p(t)=1-tq(t)$, and $\{p\in P_m:p(0)=1\}$ is exactly the image of $\{q:\mathrm{deg}q<m\}$ under $q\mapsto 1-tq$ (bijectively, the inverse being $q=(1-p)/t$, polynomial since $p(0)=1$). Minimising $\Vert b-Ax{\Vert }_2$ over ${\mathcal{K}}_m$ equals ${\min }_{p\in P_m,p(0)=1}\Vert p(A)r_0{\Vert }_2$. For the bound, $\Vert p(A)r_0{\Vert }_2\le \Vert p(A){\Vert }_2\Vert r_0{\Vert }_2$ and $p(A)=Vp(\Lambda )V^{-1}$ gives $\Vert p(A){\Vert }_2\le {\kappa }_2(V)\Vert p(\Lambda ){\Vert }_2={\kappa }_2(V){\max }_i|p({\lambda }_i)|$ since $p(\Lambda )$ is diagonal. Taking the minimum over admissible $p$ and dividing by $\Vert r_0{\Vert }_2$ yields the stated inequality. $\square$

Proposition 6 (restarted GMRES need not be monotone across cycles in the error, and may stagnate). Full GMRES residuals are monotone and reach zero by stage $\nu$; restarted GMRES$(k)$ residuals are monotone within each cycle of $k$ stages but the method can stagnate, with $\Vert r\Vert$ constant across many cycles when $k<\nu$.

Proof sketch. Within one cycle, GMRES$(k)$ is full GMRES on the current residual, so Proposition 2 gives monotone decrease across its $k$ stages. At a restart the basis $Q_k$ is discarded and Arnoldi begins afresh from the new residual $r_0^{\prime }=b-A{x}_k$, whose Krylov subspace ${\mathcal{K}}_k(A,r_0^{\prime })$ need not contain the previously useful directions. If the optimal degree-$(>k)$ residual polynomial requires roots that no degree-$k$ polynomial can supply — as when the spectrum has more than $k$ well-separated clusters, or when $A$ is non-normal so that the field of values straddles the origin — then each restarted cycle achieves the same modest reduction and the residual plateaus. The any-curve theorem (Theorem 4) shows such stagnation profiles are realisable. Stagnation is the price of the bounded storage $O(kn)$; remedies include larger $k$, preconditioning to cluster the spectrum, and augmenting the restart subspace with retained approximate eigenvectors (deflated/augmented GMRES). $\square$

Connections Master

• Krylov subspaces, Arnoldi, and Lanczos 43.07.02 is the engine GMRES runs on: the orthonormal basis $Q_m$ and the Hessenberg factorisation $A{Q}_m=Q_{m+1}{\tilde{H}}_m$ built there are exactly what turns the $n$-dimensional residual minimisation into the small $(m+1)\times m$ least-squares problem solved here. The Ritz-value residual identity $\Vert Ay-\theta y\Vert =|h_{m+1,m}||e_m^{\ast }s|$ of that unit and the GMRES residual identity $\Vert r_m\Vert =|g_{m+1}|$ are two readings of the same subdiagonal $h_{m+1,m}$ — eigenvalue accuracy there, linear-solve accuracy here — and the lucky breakdown $h_{m+1,m}=0$ means exact invariance in both, exact eigenpairs there, exact solution here.

• Least squares: normal equations, QR, and conditioning 43.04.01 supplies the inner solve: the $(m+1)\times m$ Hessenberg least-squares problem ${\min }_y\Vert \beta e_1-{\tilde{H}}_{m}y\Vert$ is solved by the incremental Givens QR of that unit, one plane rotation per stage zeroing one subdiagonal, with the last rotated component of $\beta e_1$ giving the residual norm for free. The conditioning theory there is what warns that GMRES minimises the residual, not the error: a small $\Vert b-A{x}_m\Vert$ controls the error only up to $\kappa (A)$, the very condition number that unit analyses.

• The conjugate gradient method 43.07.04 is the symmetric-positive-definite specialisation of the Krylov-solver idea developed here: where GMRES needs the full growing orthonormal basis $Q_m$ to minimise the residual for a general nonsymmetric $A$, conjugate gradients exploits the short three-term Lanczos recurrence available for $A=A^{\ast }$ to minimise the energy norm of the error with $O(n)$ storage and no stored basis. Both inherit the Krylov-subspace optimality and the polynomial-approximation convergence picture, but CG's spectrum is real and its eigenvector matrix unitary, so for CG the spectrum determines the convergence curve — the any-curve pathology of GMRES (Theorem 4) cannot occur — and the field-of-values machinery collapses back to the eigenvalue interval.

• Stationary iterative methods: Jacobi, Gauss-Seidel, SOR 43.07.01 returns as the preconditioner: a splitting $A=M-N$ that was a weak solver on its own — its spectral radius creeping to one as the grid refines — becomes the operator $M^{-1}$ that GMRES applies to cluster $\sigma (M^{-1}A)$ near one or shrink its field of values, converting the slow $1-O(h)$ convergence of even optimal SOR into the handful of stages a well-clustered GMRES needs. The clustering heuristic of Theorem 2 is precisely the design target of that preconditioning.

Historical & philosophical context Master

GMRES was introduced by Yousef Saad and Martin H. Schultz in 1986, in a paper that gave the nonsymmetric linear-system problem its first robust minimal-residual Krylov method with a guaranteed non-increasing residual and finite termination ^{[Saad, Y. & Schultz, M. H. — GMRES: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems]}. It superseded a generation of earlier nonsymmetric attempts — the unstable variants of the biconjugate gradient family and the ORTHODIR/ORTHOMIN methods — by anchoring the iterate to the numerically reliable Arnoldi orthogonalisation rather than to the short but fragile two-sided recurrences. The restart device GMRES$(k)$ appeared in the same paper as the concession to finite memory, and the stagnation it can suffer has driven the subsequent literature on deflation, augmentation, and flexible preconditioning.

The convergence theory matured more slowly and more surprisingly. The naive expectation — inherited from the symmetric conjugate-gradient case — was that clustered eigenvalues away from the origin guarantee fast convergence. Anne Greenbaum, Vlastimil Pták, and Zdeněk Strakoš showed in 1996 that this expectation is false for nonsymmetric matrices in the strongest possible sense: any non-increasing residual curve is attainable by some matrix with any prescribed spectrum, so the eigenvalues alone carry no information about the GMRES history ^{[Greenbaum, A. — Iterative Methods for Solving Linear Systems (SIAM, 1997)]}. This redirected the analysis toward the field of values, for which Howard Elman's 1982 thesis bound and later refinements give genuine rates, and toward the pseudospectra developed by Lloyd N. Trefethen and Mark Embree, whose contour-integral estimates bound the residual polynomial on a region enclosing the spectrum at a controlled distance from the origin ^{[Greenbaum, A. — Iterative Methods for Solving Linear Systems (SIAM, 1997)]}. The textbook consolidation of the algorithm, its preconditioned and restarted forms, and the polynomial-approximation convergence framework is in Saad's Iterative Methods for Sparse Linear Systems ^{[Saad, Y. — Iterative Methods for Sparse Linear Systems (2nd ed.)]}.

Bibliography Master

@article{SaadSchultz1986,
  author  = {Saad, Yousef and Schultz, Martin H.},
  title   = {{GMRES}: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems},
  journal = {SIAM Journal on Scientific and Statistical Computing},
  volume  = {7},
  number  = {3},
  year    = {1986},
  pages   = {856--869}
}

@book{Saad2003,
  author    = {Saad, Yousef},
  title     = {Iterative Methods for Sparse Linear Systems},
  edition   = {2},
  publisher = {Society for Industrial and Applied Mathematics},
  year      = {2003}
}

@book{Greenbaum1997,
  author    = {Greenbaum, Anne},
  title     = {Iterative Methods for Solving Linear Systems},
  series    = {Frontiers in Applied Mathematics},
  publisher = {SIAM},
  year      = {1997}
}

@article{GreenbaumPtakStrakos1996,
  author  = {Greenbaum, Anne and Pt{\'a}k, Vlastimil and Strako{\v{s}}, Zden{\v{e}}k},
  title   = {Any Nonincreasing Convergence Curve Is Possible for {GMRES}},
  journal = {SIAM Journal on Matrix Analysis and Applications},
  volume  = {17},
  number  = {3},
  year    = {1996},
  pages   = {465--469}
}

@book{TrefethenEmbree2005,
  author    = {Trefethen, Lloyd N. and Embree, Mark},
  title     = {Spectra and Pseudospectra: The Behavior of Nonnormal Matrices and Operators},
  publisher = {Princeton University Press},
  year      = {2005}
}

@book{TrefethenBau1997,
  author    = {Trefethen, Lloyd N. and Bau, David},
  title     = {Numerical Linear Algebra},
  publisher = {SIAM},
  address   = {Philadelphia},
  year      = {1997}
}

@phdthesis{Elman1982,
  author = {Elman, Howard C.},
  title  = {Iterative Methods for Large, Sparse, Nonsymmetric Systems of Linear Equations},
  school = {Yale University},
  year   = {1982}
}