Right preconditioned GMRES for arbitrary singular systems

Brown and Walker (1997) showed that GMRES determines a least squares solution of $A x = b$ where $A \in {\bf R}^{n \times n}$ without breakdown for arbitrary $b, x_0 \in {\bf R}^n$ if and only if $A$ is range-symmetric, i.e. ${\cal R} (A^{\rm T}) = {\cal R} (A)$, where $A$ may be singular and $b$ may not be in the range space ${\cal R} A)$ of $A$. In this paper, we propose applying GMRES to $A C A^{\rm T} z = b$, where $C \in {\bf R}^{n \times n}$ is symmetric positive definite. This determines a least squares solution $x = CA^{\rm T} z$ of $A x = b$ without breakdown for arbitrary (singular) matrix $A \in {\bf R}^{n \times n}$ and $b \in {\bf R}^n$. To make the method numerically stable, we propose using the pseudoinverse with an appropriate threshold parameter to suppress the influence of tiny singular values when solving the severely ill-conditioned Hessenberg systems which arise in the Arnoldi process of GMRES when solving inconsistent range-symmetric systems. Numerical experiments show that the method taking $C$ to be the identity matrix and the inverse matrix of the diagonal matrix whose diagonal elements are the diagonal of $A A^{\rm T}$ gives a least squares solution even when $A$ is not range-symmetric, including the case when ${\rm index}(A) >1$.