Iterated Gauss-Seidel GMRES

The GMRES algorithm of Saad and Schultz (1986) is an iterative method for approximately solving linear systems $A{\bf x}={\bf b}$, with initial guess ${\bf x}_0$ and residual ${\bf r}_0 = {\bf b} - A{\bf x}_0$. The algorithm employs the Arnoldi process to generate the Krylov basis vectors (the columns of $V_k$). It is well known that this process can be viewed as a $QR$ factorization of the matrix $B_k = [\: {\bf r}_0, AV_k\:]$ at each iteration. Despite an ${\cal O}(\epsilon)\kappa(B_k)$ loss of orthogonality, for unit roundoff $\epsilon$ and condition number $\kappa$, the modified Gram-Schmidt formulation was shown to be backward stable in the seminal paper by Paige et al.~(2006). We present an iterated Gauss-Seidel formulation of the GMRES algorithm (IGS-GMRES) based on the ideas of Ruhe (1983) and \'{S}wirydowicz et al. (2020). IGS-GMRES maintains orthogonality to the level ${\cal O}(\epsilon)\kappa(B_k)$ or ${\cal O}(\epsilon)$, depending on the choice of one or two iterations; for two Gauss-Seidel iterations, the computed Krylov basis vectors remain linearly independent and the smallest singular value of $V_k$ remains close to one. The resulting GMRES method is thus backward stable. We show that IGS-GMRES can be implemented with only a single synchronization point per iteration, making it relevant to large-scale parallel computing environments. We also demonstrate that, unlike MGS-GMRES, in IGS-GMRES the relative Arnoldi residual corresponding to the computed approximate solution no longer stagnates above machine precision even for highly non-normal systems.