Documentation

gmres

Generalized minimum residual method (with restarts)

collapse all in page

Syntax

x = gmres(A,b)
gmres(A,b,restart)
gmres(A,b,restart,tol)
gmres(A,b,restart,tol,maxit)
gmres(A,b,restart,tol,maxit,M)
gmres(A,b,restart,tol,maxit,M1,M2)
gmres(A,b,restart,tol,maxit,M1,M2,x0)
[x,flag] = gmres(A,b,...)
[x,flag,relres] = gmres(A,b,...)
[x,flag,relres,iter] = gmres(A,b,...)
[x,flag,relres,iter,resvec] = gmres(A,b,...)

Description

x = gmres(A,b)attempts to solve the system of linear equationsA*x = bforx. Then-by-ncoefficient matrixAmust be square and should be large and sparse. The column vectorbmust have lengthn.Acan be a function handle,afun, such thatafun(x)returnsA*x. For this syntax,gmresdoes not restart; the maximum number of iterations ismin(n,10).

Parameterizing Functionsexplains how to provide additional parameters to the functionafun, as well as the preconditioner functionmfundescribed below, if necessary.

Ifgmresconverges, a message to that effect is displayed. Ifgmresfails to converge after the maximum number of iterations or halts for any reason, a warning message is printed displaying the relative residualnorm(b-A*x)/norm(b)and the iteration number at which the method stopped or failed.

gmres(A,b,restart)restarts the method everyrestartinner iterations. The maximum number of outer iterations ismin(n/restart,10). The maximum number of total iterations isrestart*min(n/restart,10). Ifrestartisnor[], thengmresdoes not restart and the maximum number of total iterations ismin(n,10).

gmres(A,b,restart,tol)specifies the tolerance of the method. Iftolis[], thengmresuses the default,1e-6.

gmres(A,b,restart,tol,maxit)specifies the maximum number of outer iterations, i.e., the total number of iterations does not exceedrestart*maxit. Ifmaxitis[]thengmresuses the default,min(n/restart,10). Ifrestartisnor[], then the maximum number of total iterations ismaxit(instead ofrestart*maxit).

gmres(A,b,restart,tol,maxit,M)andgmres(A,b,restart,tol,maxit,M1,M2)use preconditionerMorM = M1*M2and effectively solve the systeminv(M)*A*x = inv(M)*bforx. IfMis[]thengmresapplies no preconditioner.Mcan be a function handlemfunsuch thatmfun(x)returnsM\x.

gmres(A,b,restart,tol,maxit,M1,M2,x0)specifies the first initial guess. Ifx0is[], thengmresuses the default, an all-zero vector.

[x,flag] = gmres(A,b,...)also returns a convergence flag:

flag = 0

gmresconverged to the desired tolerancetolwithinmaxitouter iterations.

flag = 1

gmresiteratedmaxittimes but did not converge.

flag = 2

PreconditionerMwas ill-conditioned.

flag = 3

gmresstagnated. (Two consecutive iterates were the same.)

Wheneverflagis not0, the solutionxreturned is that with minimal norm residual computed over all the iterations. No messages are displayed if theflagoutput is specified.

[x,flag,relres] = gmres(A,b,...)also returns the relative residualnorm(b-A*x)/norm(b). Ifflagis0,relres <= tol. The third output,relres, is the relative residual of the preconditioned system.

[x,flag,relres,iter] = gmres(A,b,...)also returns both the outer and inner iteration numbers at whichxwas computed, where0 <= iter(1) <= maxitand0 <= iter(2) <= restart.

[x,flag,relres,iter,resvec] = gmres(A,b,...)also returns a vector of the residual norms at each inner iteration. These are the residual norms for the preconditioned system.

Examples

使用gmres with a Matrix Input

A = gallery('wilk',21); b = sum(A,2); tol = 1e-12; maxit = 15; M1 = diag([10:-1:1 1 1:10]); x = gmres(A,b,10,tol,maxit,M1);

显示以下信息:

gmr(10)聚集在外层迭代2(内eration 9) to a solution with relative residual 3.3e-013

使用gmres with a Function Handle

This example replaces the matrixAin the previous example with a handle to a matrix-vector product functionafun, and the preconditionerM1with a handle to a backsolve functionmfun. The example is contained in a functionrun_gmresthat

  • Callsgmreswith the function handle@afunas its first argument.

  • Containsafunandmfunas nested functions, so that all variables inrun_gmresare available toafunandmfun.

The following shows the code forrun_gmres:

function x1 = run_gmres n = 21; b = afun(ones(n,1)); tol = 1e-12; maxit = 15; x1 = gmres(@afun,b,10,tol,maxit,@mfun); function y = afun(x) y = [0; x(1:n-1)] + ... [((n-1)/2:-1:0)'; (1:(n-1)/2)'].*x + ... [x(2:n); 0]; end function y = mfun(r) y = r ./ [((n-1)/2:-1:1)'; 1; (1:(n-1)/2)']; end end

When you enter

x1 = run_gmres;

MATLAB®software displays the message

gmr(10)聚集在外层迭代2(内eration 10) to a solution with relative residual 1.1e-013.

使用a Preconditioner without Restart

This example demonstrates the use of a preconditioner without restartinggmres.

Loadwest0479, a real 479-by-479 nonsymmetric sparse matrix.

loadwest0479; A = west0479;

Set the tolerance and maximum number of iterations.

tol = 1e-12; maxit = 20;

Definebso that the true solution is a vector of all ones.

b = full(sum(A,2)); [x0,fl0,rr0,it0,rv0] = gmres(A,b,[],tol,maxit);

fl0is 1 becausegmres不收敛到所请求的宽容1e-12within the requested 20 iterations. The best approximate solution thatgmresreturns is the last one (as indicated byit0(2) = 20). MATLAB stores the residual history inrv0.

Plot the behavior ofgmres.

semilogy(0:maxit,rv0/norm(b),'-o'); xlabel('Iteration number'); ylabel('Relative residual');

The plot shows that the solution converges slowly. A preconditioner may improve the outcome.

Useiluto form the preconditioner, sinceAis nonsymmetric.

[L,U] = ilu(A,struct('type','ilutp','droptol',1e-5));
Error using ilu There is a pivot equal to zero. Consider decreasing the drop tolerance or consider using the 'udiag' option.

Note MATLAB cannot construct the incomplete LU as it would result in a singular factor, which is useless as a preconditioner.

As indicated by the error message, try again with a reduced drop tolerance.

[L,U] = ilu(A,struct('type','ilutp','droptol',1e-6)); [x1,fl1,rr1,it1,rv1] = gmres(A,b,[],tol,maxit,L,U);

fl1is 0 becausegmresdrives the relative residual to9.5436e-14(the value ofrr1). The relative residual is less than the prescribed tolerance of1e-12at the sixth iteration (the value ofit1(2)) when preconditioned by the incomplete LU factorization with a drop tolerance of1e-6. The output,rv1(1), isnorm(M\b), whereM = L*U. The output,rv1(7), isnorm(U\(L\(b-A*x1))).

Follow the progress ofgmresby plotting the relative residuals at each iteration starting from the initial estimate (iterate number 0).

semilogy(0:it1(2),rv1/norm(b),'-o'); xlabel('Iteration number'); ylabel('Relative residual');

使用a Preconditioner with Restart

This example demonstrates the use of a preconditioner with restartedgmres.

Loadwest0479, a real 479-by-479 nonsymmetric sparse matrix.

loadwest0479; A = west0479;

Definebso that the true solution is a vector of all ones.

b = full(sum(A,2));

Construct an incomplete LU preconditioner as in the previous example.

[L,U] = ilu(A,struct('type','ilutp','droptol',1e-6));

The benefit to using restartedgmresis to limit the amount of memory required to execute the method. Without restart,gmresrequiresmaxitvectors of storage to keep the basis of the Krylov subspace. Also,gmresmust orthogonalize against all of the previous vectors at each step. Restarting limits the amount of workspace used and the amount of work done per outer iteration. Note that even though preconditionedgmresconverged in six iterations above, the algorithm allowed for as many as twenty basis vectors and therefore, allocated all of that space up front.

Executegmres(3),gmres(4), andgmres(5)

tol = 1e-12; maxit = 20; re3 = 3; [x3,fl3,rr3,it3,rv3] = gmres(A,b,re3,tol,maxit,L,U); re4 = 4; [x4,fl4,rr4,it4,rv4] = gmres(A,b,re4,tol,maxit,L,U); re5 = 5; [x5,fl5,rr5,it5,rv5] = gmres(A,b,re5,tol,maxit,L,U);

fl3,fl4, andfl5are all 0 because in each case restartedgmresdrives the relative residual to less than the prescribed tolerance of1e-12.

The following plots show the convergence histories of each restartedgmresmethod.gmres(3)converges at outer iteration 5, inner iteration 3 (it3 = [5, 3]) which would be the same as outer iteration 6, inner iteration 0, hence the marking of 6 on the final tick mark.

figure semilogy(1:1/3:6,rv3/norm(b),'-o'); h1 = gca; h1.XTick = [1:1/3:6]; h1.XTickLabel = ['1';' ';' ';'2';' ';' ';'3';' ';' ';'4';' ';' ';'5';' ';' ';'6';]; title('gmres(3)') xlabel('Iteration number'); ylabel('Relative residual');

figure semilogy(1:1/4:3,rv4/norm(b),'-o'); h2 = gca; h2.XTick = [1:1/4:3]; h2.XTickLabel = ['1';' ';' ';' ';'2';' ';' ';' ';'3']; title('gmres(4)') xlabel('Iteration number'); ylabel('Relative residual');

figure semilogy(1:1/5:2.8,rv5/norm(b),'-o'); h3 = gca; h3.XTick = [1:1/5:2.8]; h3.XTickLabel = ['1';' ';' ';' ';' ';'2';' ';' ';' ';' ']; title('gmres(5)') xlabel('Iteration number'); ylabel('Relative residual');

References

Barrett, R., M. Berry, T. F. Chan, et al.,Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods, SIAM, Philadelphia, 1994.

Saad, Youcef and Martin H. Schultz, “GMRES: A generalized minimal residual algorithm for solving nonsymmetric linear systems,”SIAM J. Sci. Stat. Comput., July 1986, Vol. 7, No. 3, pp. 856-869.

See Also

|||||||||

Topics

Introduced before R2006a

Was this topic helpful?
MATLAB Documentation
Other Documentation
金宝app

Try MATLAB, Simulink, and Other Products

Get trial now