GMRES solves the problem $Ax = b$ approximately for $x$ where $A$ is a general, linear operator and $b$ the right-hand side vector. The method is optimal in the sense that it selects the solution with minimal residual from a Krylov subspace, but the price of optimality is increasing storage and computation effort per iteration. Restarts are necessary to fix these costs.
gmres(A, b; kwargs...) -> x, [history]
gmres!, but allocates a solution vector
x initialized with zeros.
gmres!(x, A, b; kwargs...) -> x, [history]
Solves the problem $Ax = b$ with restarted GMRES.
x: Initial guess, will be updated in-place;
A: linear operator;
b: right-hand side.
iszero(x)so that one matrix-vector product can be saved when computing the initial residual vector;
tol: relative tolerance;
restart::Int = min(20, size(A, 2)): restarts GMRES after specified number of iterations;
maxiter::Int = size(A, 2): maximum number of inner iterations of GMRES;
Pl: left preconditioner;
Pr: right preconditioner;
log::Bool: keep track of the residual norm in each iteration;
verbose::Bool: print convergence information during the iterations.
x: approximate solution.
x: approximate solution;
history: convergence history.
The implementation pre-allocates a matrix $V$ of size
restart whose columns form an orthonormal basis for the Krylov subspace. This allows BLAS2 operations when updating the solution vector $x$. The Hessenberg matrix is also pre-allocated.
Modified Gram-Schmidt is used to orthogonalize the columns of $V$.
The computation of the residual norm is implemented in a non-standard way, namely keeping track of a vector $\gamma$ in the null-space of $H_k^*$, which is the adjoint of the $(k + 1) \times k$ Hessenberg matrix $H_k$ at the $k$th iteration. Only when $x$ needs to be updated is the Hessenberg matrix mutated with Givens rotations.
GMRES can be used as an iterator. This makes it possible to access the Hessenberg matrix and Krylov basis vectors during the iterations.