《gmres源程序 matlab.doc》由会员分享，可在线阅读，更多相关《gmres源程序 matlab.doc（14页珍藏版）》请在金锄头文库上搜索。
1、function x,flag,relres,iter,resvec = gmres(A,b,restart,tol,maxit,M1,M2,x,varargin)%GMRES Generalized Minimum Residual Method.% X = GMRES(A,B) attempts to solve the system of linear equations A*X = B% for X. The N-by-N coefficient matrix A must be square and the right% hand side column vector B must
2、have length N. This uses the unrestarted% method with MIN(N,10) total iterations.% X = GMRES(AFUN,B) accepts a function handle AFUN instead of the matrix% A. AFUN(X) accepts a vector input X and returns the matrix-vector% product A*X. In all of the following syntaxes, you can replace A by% AFUN.% X
3、= GMRES(A,B,RESTART) restarts the method every RESTART iterations.% If RESTART is N or then GMRES uses the unrestarted method as above.% X = GMRES(A,B,RESTART,TOL) specifies the tolerance of the method. If% TOL is then GMRES uses the default, 1e-6.% X = GMRES(A,B,RESTART,TOL,MAXIT) specifies the max
4、imum number of outer% iterations. Note: the total number of iterations is RESTART*MAXIT. If% MAXIT is then GMRES uses the default, MIN(N/RESTART,10). If RESTART% is N or then the total number of iterations is MAXIT.% X = GMRES(A,B,RESTART,TOL,MAXIT,M) and% X = GMRES(A,B,RESTART,TOL,MAXIT,M1,M2) use
5、preconditioner M or M=M1*M2% and effectively solve the system inv(M)*A*X = inv(M)*B for X. If M is% then a preconditioner is not applied. M may be a function handle% returning MX.% X = GMRES(A,B,RESTART,TOL,MAXIT,M1,M2,X0) specifies the first initial% guess. If X0 is then GMRES uses the default, an
6、all zero vector.% X,FLAG = GMRES(A,B,.) also returns a convergence FLAG:% 0 GMRES converged to the desired tolerance TOL within MAXIT iterations.% 1 GMRES iterated MAXIT times but did not converge.% 2 preconditioner M was ill-conditioned.% 3 GMRES stagnated (two consecutive iterates were the same).%
7、 X,FLAG,RELRES = GMRES(A,B,.) also returns the relative residual% NORM(B-A*X)/NORM(B). If FLAG is 0, then RELRES = TOL. Note with% preconditioners M1,M2, the residual is NORM(M2(M1(B-A*X).% X,FLAG,RELRES,ITER = GMRES(A,B,.) also returns both the outer and% inner iteration numbers at which X was comp
8、uted: 0 = ITER(1) = MAXIT% and 0 = ITER(2) = RESTART.% X,FLAG,RELRES,ITER,RESVEC = GMRES(A,B,.) also returns a vector of% the residual norms at each inner iteration, including NORM(B-A*X0).% Note with preconditioners M1,M2, the residual is NORM(M2(M1(B-A*X).% Example:% n = 21; A = gallery(wilk,n); b
9、 = sum(A,2);% tol = 1e-12; maxit = 15; M = diag(10:-1:1 1 1:10);% x = gmres(A,b,10,tol,maxit,M);% Or, use this matrix-vector product function% %-% function y = afun(x,n)% y = 0; x(1:n-1) + (n-1)/2:-1:0); (1:(n-1)/2).*x+x(2:n); 0;% %-% and this preconditioner backsolve function% %-% function y = mfun
10、(r,n)% y = r ./ (n-1)/2:-1:1); 1; (1:(n-1)/2);% %-% as inputs to GMRES:% x1 = gmres(x)afun(x,n),b,10,tol,maxit,(x)mfun(x,n);% Class support for inputs A,B,M1,M2,X0 and the output of AFUN:% float: double% See also BICG, BICGSTAB, BICGSTABL, CGS, LSQR, MINRES, PCG, QMR, SYMMLQ,% TFQMR, ILU, FUNCTION_H
11、ANDLE.% References% H.F. Walker, Implementation of the GMRES Method Using Householder% Transformations, SIAM J. Sci. Comp. Vol 9. No 1. January 1988.% Copyright 1984-2010 The MathWorks, Inc.% $Revision: 1.21.4.13 $ $Date: 2010/02/25 08:11:48 $if (nargin 2) error(MATLAB:gmres:NumInputs,Not enough inp
12、ut arguments.);end% Determine whether A is a matrix or a function.atype,afun,afcnstr = iterchk(A);if strcmp(atype,matrix) % Check matrix and right hand side vector inputs have appropriate sizes m,n = size(A); if (m = n) error(MATLAB:gmres:SquareMatrix,Matrix must be square.); end if isequal(size(b),
13、m,1) error(MATLAB:gmres:VectorSize, %s %d %s, . Right hand side must be a column vector of length, m,. to match the coefficient matrix.); endelse m = size(b,1); n = m; if iscolumn(b) error(MATLAB:gmres:Vector,Right hand side must be a column vector.); endend% Assign default values to unspecified parametersif (nargin 3) | isempty(restart) | (restart = n) restarted = false;else restarted = true;endif (nargin 4) | isempty(tol) tol = 1e-6;endwarned = 0;if tol eps warning(MATLAB:gmres:tooSmallTolerance, . strcat(Input tol is smaller than eps a
