Flexible GMRES with variable right preconditioning

Let $n,m\in \mathbb{N}$ with $n,m\ge 1$ (The natural numbers $\mathbb{N}$ (von Neumann)), let $A\in M_n(\mathbb{C})$ (The vector space $M_{m\times n}(F):=F^{m\times n}$ of $m$ by $n$ matrices over a field, with entrywise operations), and let $b,x_0\in {\mathbb{C}}^n$. The complex scalars form a field ($\mathbb{C}=\mathbb{R}[x]\mathrm{/}(x^2+1)$ is a field, every element is uniquely $a+bi$, and every nonzero element has inverse $(a-bi)\mathrm{/}(a^2+b^2)$). As additional input, take maps ${\mathcal{P}}_1,\dots ,{\mathcal{P}}_m:{\mathbb{C}}^n\to {\mathbb{C}}^n$. The maps may vary with $j$; no linearity or fixedness is assumed. In the usual matrix case, ${\mathcal{P}}_j(v)=M_j^{-1}v$ for an invertible $M_j$ (Invertible matrices and the general linear group ${\mathrm{GL}}_n(F)$). Equip each coordinate space ${\mathbb{C}}^q$ used below with the standard Hermitian inner product $⟨z,w⟩={\sum }_{a<q}{z}_a\overline{w_a}$ and its induced Euclidean norm $\mathrm{\parallel }z{\mathrm{\parallel }}_2=\sqrt{⟨z,z⟩}$ (The standard formulas $⟨x,y⟩={\sum }_{k<n}{x}_k{y}_k$ on ${\mathbb{R}}^n$ and ${\sum }_{k<n}{x}_k\overline{y_k}$ on ${\mathbb{C}}^n$ are inner products, The norm $\parallel v\parallel =\sqrt{⟨v,v⟩}$ induced by a real or complex inner product). All arithmetic and orthogonalization are exact. Iteration and auxiliary Hessenberg labels below start at $1$; coordinates of $A$ and ambient vectors retain the library's zero-based convention.

The algorithm uses and stores each actual $z_j$ and makes no fixed-matrix assumption.

Set $r_0=b-A{x}_0$ and $\beta =\mathrm{\parallel }{r}_0{\mathrm{\parallel }}_2$. If $\beta =0$, return initial-solution($x_0$). Otherwise set $v_1=r_0\mathrm{/}\beta$. For $j=1,\dots ,m$, perform:

1. Obtain $z_j={\mathcal{P}}_j(v_j)$ and set $w:=A{z}_j$.
2. For $i=1,\dots ,j$ in order, set $h_{ij}:=⟨w,v_i⟩,w:=w-h_{ij}{v}_i\mathrm{.}$
3. Set $h_{j+1,j}:=\mathrm{\parallel }w{\mathrm{\parallel }}_2$. If $h_{j+1,j}=0$, return breakdown($j,x_{j-1}$) with the last fully formed iterate (for $j=1$ it is $x_0$); do not normalize or claim a new iterate at this step.
4. If $h_{j+1,j}>0$, set $v_{j+1}:=w\mathrm{/}{h}_{j+1,j}$. Let $e_0\in {\mathbb{C}}^{j+1}$ be the standard unit vector at coordinate $0$ (coordinate $0$ is $1$ and all other coordinates are $0$; The standard list $e:n\to F^n$ with $e_i(i)=1_F$ and $e_i(j)=0_F$ for $j\mathrm{\ne }i$ is an ordered basis of $F^n$; hence ${\dim }_F{F}^n=n$, and $F^0$ is the zero space with basis $\mathrm{\varnothing }$ and dimension $0$). The iteration labels $h_{ik}$ remain one-based, while matrix slots are zero-based: put $Z_j\in M_{n\times j}(\mathbb{C})$ with $(Z_j{)}_{a,k-1}=(z_k{)}_a$ for $0\le a<n$, $1\le k\le j$, and let ${\bar{H}}_j\in M_{(j+1)\times j}(\mathbb{C})$ satisfy $({\bar{H}}_j{)}_{i-1,k-1}=h_{ik}$ for $1\le k\le j$, $1\le i\le k+1$, and $({\bar{H}}_j{)}_{a,b}=0$ when $a>b+1$. Choose any $y_j\in {\mathrm{argmin}}_{y\in {\mathbb{C}}^j}\mathrm{\parallel }\beta e_0-{\bar{H}}_{j}y{\mathrm{\parallel }}_2$ and set $x_j:=x_0+Z_j{y}_j$.

If all $m$ steps complete without breakdown, return budget-complete($m,x_m,Z_m$). The finite-dimensional least-squares minimizer exists because each coordinate space ${\mathbb{C}}^q$ is finite-dimensional by The standard list $e:n\to F^n$ with $e_i(i)=1_F$ and $e_i(j)=0_F$ for $j\mathrm{\ne }i$ is an ordered basis of $F^n$; hence ${\dim }_F{F}^n=n$, and $F^0$ is the zero space with basis $\mathrm{\varnothing }$ and dimension $0$, and the least-squares theorem For a linear map $T:V\to W$ between finite-dimensional inner-product spaces, $x$ minimises $\parallel Tx-b\parallel$ if and only if $T^{\ast }(Tx-b)=0$, equivalently $T^{\ast }Tx=T^{\ast }b$; minimisers exist and any two differ by an element of $\ker T$ applies; the output relation permits any minimizer. The algorithm has no convergence test: budget-complete and breakdown assert neither convergence nor that a breakdown iterate solves $Ax=b$. The residual-minimization property over the actual span of the $z_j$ is proved separately in Flexible GMRES minimizes over the actual preconditioned-vector span. Every division is by a positive real norm: the initial $\beta =0$ case returns before normalizing $r_0$, and each later $h_{j+1,j}=0$ case returns before normalizing $w$.