Minimal residual method

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.{{cite thesis |first=M. Naoufal |last=Nifa |title=Efficient solvers for constrained optimization in parameter identification problems|date=24 November 2017 | type=Doctoral Thesis|url=https://www.theses.fr/2017SACLC066.pdf|pages=51–52|publisher=Université Paris Saclay (COmUE) }}

The MINRES method iteratively calculates an approximate solution of a linear system of equations of the form

$$Ax\; =\; b,$$

where $A\backslash in\backslash mathbb\{R\}^\{n\backslash times\; n\}$ is a symmetric matrix and $b\backslash in\backslash mathbb\{R\}^n$ a vector.

For this, the norm of the residual $r(x)\; :=\; b\; -\; Ax$ in a $k$-dimensional Krylov subspace

$$V\_k\; =\; x\_0\; +\; \backslash operatorname\{span\}\backslash \{r\_0,\; Ar\_0\backslash ldots,A^\{k-1\}r\_0\backslash \}$$

is minimized. Here $x\_0\backslash in\backslash mathbb\{R\}^n$ is an initial value (often zero) and $r\_0\; :=\; r(x\_0)$.

More precisely, we define the approximate solutions $x\_k$ through

$$x\_k\; :=\; \backslash mathrm\{argmin\}\_\{x\backslash in\; V\_k\}\; \backslash |r(x)\backslash |,$$

where $\backslash |\backslash cdot\backslash |$ is the standard Euclidean norm on $\backslash mathbb\{R\}^n$.

Because of the symmetry of $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.

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 $p\_k$) can be orthogonalized, whereas in CR their images (below labeled with $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 $x\_0\backslash in\backslash mathbb\{R\}^n$ arbitrary and compute

$$\backslash begin\{align\}$$

r_0 &= b - A x_0 \\

p_0 &= r_0 \\

s_0 &= A p_0

\end{align}

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

{{unordered list

| Compute $x\_k,r\_k$ through

$$\backslash alpha\_\{k-1\}\; =\; \backslash frac\{\; \backslash langle\; r\_\{k-1\},\; s\_\{k-1\}\; \backslash rangle\}\{\; \backslash langle\; s\_\{k-1\},\; s\_\{k-1\}\; \backslash rangle\; \}$$

$$x\_k\; =\; x\_\{k-1\}\; +\; \backslash alpha\_\{k-1\}\; p\_\{k-1\}$$

$$r\_k\; =\; r\_\{k-1\}\; -\; \backslash alpha\_\{k-1\}\; s\_\{k-1\}$$

if $\backslash |r\_k\backslash |$ is smaller than a specified tolerance, the algorithm is interrupted with the approximate solution $x\_k$. Otherwise, a new descent direction $p\_k$ is calculated through

$$p\_k\; \backslash leftarrow\; s\_\{k-1\}$$

$$s\_k\; \backslash leftarrow\; As\_\{k-1\}$$

| for $l=1,2$ (the step $l=2$ is not carried out in the first iteration step) calculate:

$$\backslash beta\_\{k,l\}\; =\; \backslash frac\{\backslash langle\; s\_k,\; s\_\{k-l\}\; \backslash rangle\}\{\; \backslash langle\; s\_\{k-l\},\; s\_\{k-l\}\; \backslash rangle\}$$

$$p\_k\; \backslash leftarrow\; p\_k\; -\; \backslash beta\_\{k,l\}\; p\_\{k-l\}$$

$$s\_k\; \backslash leftarrow\; s\_k\; -\; \backslash beta\_\{k,l\}\; s\_\{k-l\}$$

}}

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.{{Cite book |title=Numerical Methods for Two-phase Incompressible Flows |publisher=Springer | author=Sven Gross, Arnold Reusken |date=6 May 2011 |location=section 5.2 |isbn=978-3-642-19685-0}} 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:

$$\backslash |r\_k\backslash |\; \backslash le\; 2\backslash left(\backslash frac\{\backslash sqrt\{\backslash kappa(A)\}-1\}\{\backslash sqrt\{\backslash kappa(A)\}+1\}\backslash right)^k\backslash |r\_\{0\}\backslash |,$$

where $\backslash kappa(A)$ is the condition number of matrix $A$. Because $A$ is normal, we have

$$\backslash kappa(A)\; =\; \backslash frac\{\backslash left|\backslash lambda\_\backslash text\{max\}(A)\backslash right|\}\{\backslash left|\backslash lambda\_\backslash text\{min\}(A)\backslash right|\},$$

where $\backslash lambda\_\backslash text\{max\}(A)$ and $\backslash lambda\_\backslash text\{min\}(A)$ are maximal and minimal eigenvalues of $A$, respectively.

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

• [https://mathworld.wolfram.com/MinimalResidualMethod.html Minimal Residual Method], Wolfram MathWorld, Jul 26, 2022.

Category:Numerical linear algebra