Minimal residual method

The Minimal Residual Method or MINRES is a Krylov subspace method for the iterative solution of symmetric linear equation systems. It was proposed by mathematicians Christopher Conway Paige and Michael Alan Saunders in 1975.^[1]

A comparison of the norm of error and residual in the CG method (blue) and the MINRES method (green). The matrix used comes from a 2D boundary-value problem.

In contrast to the popular CG method, the MINRES method does not assume that the matrix is positive definite, only the symmetry of the matrix is mandatory.

GMRES vs. MINRES

The GMRES method is essentially a generalization of MINRES for arbitrary matrices. Both minimize the 2-norm of the residual and do the same calculations in exact arithmetic when the matrix is symmetric. MINRES is a short-recurrence method with a constant memory requirement, whereas GMRES requires storing the whole Krylov space, so its memory requirement is roughly proportional to the number of iterations. On the other hand, GMRES tends to suffer less from loss of orthogonality.^[1][2]

Properties of the MINRES method

The MINRES method iteratively calculates an approximate solution of a linear system of equations of the form $${\displaystyle Ax=b,}$$ where ${\displaystyle A\in \mathbb {R} ^{n\times n}}$ is a symmetric matrix and ${\displaystyle b\in \mathbb {R} ^{n}}$ a vector.

For this, the norm of the residual ${\displaystyle r(x):=b-Ax}$ in a ${\displaystyle k}$-dimensional Krylov subspace $${\displaystyle V_{k}=x_{0}+\operatorname {span} \{r_{0},Ar_{0}\ldots ,A^{k-1}r_{0}\}}$$ is minimized. Here ${\displaystyle x_{0}\in \mathbb {R} ^{n}}$ is an initial value (often zero) and ${\displaystyle r_{0}:=r(x_{0})}$.

More precisely, we define the approximate solutions ${\displaystyle x_{k}}$ through $${\displaystyle x_{k}:=\mathrm {argmin} _{x\in V_{k}}\|r(x)\|,}$$ where ${\displaystyle \|\cdot \|}$ is the standard Euclidean norm on ${\displaystyle \mathbb {R} ^{n}}$.

Because of the symmetry of ${\displaystyle A}$, unlike in the GMRES method, it is possible to carry out this minimization process recursively, storing only two previous steps (short recurrence). This saves memory.

MINRES algorithm

Note: The MINRES method is more complicated than the algebraically equivalent Conjugate Residual method. The Conjugate Residual (CR) method was therefore produced below as a substitute. It differs from MINRES in that in MINRES, the columns of a basis of the Krylov space (denoted below by ${\displaystyle p_{k}}$) can be orthogonalized, whereas in CR their images (below labeled with ${\displaystyle s_{k}}$) can be orthogonalized via the Lanczos recursion. There are more efficient and preconditioned variants with fewer AXPYs. Compare with the article.

First you choose ${\displaystyle x_{0}\in \mathbb {R} ^{n}}$ arbitrary and compute $${\displaystyle {\begin{aligned}r_{0}&=b-Ax_{0}\\p_{0}&=r_{0}\\s_{0}&=Ap_{0}\end{aligned}}}$$

Then we iterate for ${\displaystyle k=1,2,\dots }$ in the following steps:

• Compute ${\displaystyle x_{k},r_{k}}$ through $${\displaystyle \alpha _{k-1}={\frac {\langle r_{k-1},s_{k-1}\rangle }{\langle s_{k-1},s_{k-1}\rangle }}}$$ $${\displaystyle x_{k}=x_{k-1}+\alpha _{k-1}p_{k-1}}$$ $${\displaystyle r_{k}=r_{k-1}-\alpha _{k-1}s_{k-1}}$$ if ${\displaystyle \|r_{k}\|}$ is smaller than a specified tolerance, the algorithm is interrupted with the approximate solution ${\displaystyle x_{k}}$. Otherwise, a new descent direction ${\displaystyle p_{k}}$ is calculated through $${\displaystyle p_{k}\leftarrow s_{k-1}}$$ $${\displaystyle s_{k}\leftarrow As_{k-1}}$$
• for ${\displaystyle l=1,2}$ (the step ${\displaystyle l=2}$ is not carried out in the first iteration step) calculate: $${\displaystyle \beta _{k,l}={\frac {\langle s_{k},s_{k-l}\rangle }{\langle s_{k-l},s_{k-l}\rangle }}}$$ $${\displaystyle p_{k}\leftarrow p_{k}-\beta _{k,l}p_{k-l}}$$ $${\displaystyle s_{k}\leftarrow s_{k}-\beta _{k,l}s_{k-l}}$$

Convergence rate of the MINRES method

In the case of positive definite matrices, the convergence rate of the MINRES method can be estimated in a way similar to that of the CG method.^[3] In contrast to the CG method, however, the estimation does not apply to the errors of the iterates, but to the residual. The following applies:

$${\displaystyle \|r_{k}\|\leq 2\left({\frac {{\sqrt {\kappa (A)}}-1}{{\sqrt {\kappa (A)}}+1}}\right)^{k}\|r_{0}\|,}$$

where ${\displaystyle \kappa (A)}$ is the condition number of matrix ${\displaystyle A}$. Because ${\displaystyle A}$ is normal, we have $${\displaystyle \kappa (A)={\frac {\left|\lambda _{\text{max}}(A)\right|}{\left|\lambda _{\text{min}}(A)\right|}},}$$ where ${\displaystyle \lambda _{\text{max}}(A)}$ and ${\displaystyle \lambda _{\text{min}}(A)}$ are maximal and minimal eigenvalues of ${\displaystyle A}$, respectively.

Implementation in GNU Octave / MATLAB

function [x, r] = minres(A, b, x0, maxit, tol)
  x = x0;
  r = b - A * x0;
  p0 = r;
  s0 = A * p0;
  p1 = p0;
  s1 = s0;
  for iter = 1:maxit
    p2 = p1; p1 = p0;
    s2 = s1; s1 = s0;
    alpha = r'*s1 / (s1'*s1);
    x = x + alpha * p1;
    r = r - alpha * s1;
    if (r'*r < tol^2)
      break
    end
    p0 = s1;
    s0 = A * s1;
    beta1 = s0'*s1 / (s1'*s1);
    p0 = p0 - beta1 * p1;
    s0 = s0 - beta1 * s1;
    if iter > 1
      beta2 = s0'*s2 / (s2'*s2);
      p0 = p0 - beta2 * p2;
      s0 = s0 - beta2 * s2;
    end
  end
end

References

1. Christopher C. Paige, Michael A. Saunders (1975). "Solution of sparse indefinite systems of linear equations". SIAM Journal on Numerical Analysis. 12 (4): 617–629. doi:10.1137/0712047.
2. Nifa, M. Naoufal. "Efficient solvers for constrained optimization in parameter identification problems" (PDF) (Doctoral Thesis). pp. 51–52.
3. Sven Gross, Arnold Reusken. Numerical Methods for Two-phase Incompressible Flows. section 5.2: Springer. ISBN 978-3-642-19685-0.{{cite book}}: CS1 maint: location (link)

External links

• Minimal Residual Method, Wolfram MathWorld, Jul 26, 2022.