gmres源程序matlabWord格式文档下载.docx
- 文档编号:18591875
- 上传时间:2022-12-28
- 格式:DOCX
- 页数:19
- 大小:19.25KB
gmres源程序matlabWord格式文档下载.docx
《gmres源程序matlabWord格式文档下载.docx》由会员分享,可在线阅读,更多相关《gmres源程序matlabWord格式文档下载.docx(19页珍藏版)》请在冰豆网上搜索。
%[X,FLAG,RELRES]=GMRES(A,B,...)alsoreturnstherelativeresidual
%NORM(B-A*X)/NORM(B).IfFLAGis0,thenRELRES<
=TOL.Notewith
%preconditionersM1,M2,theresidualisNORM(M2\(M1\(B-A*X))).
%[X,FLAG,RELRES,ITER]=GMRES(A,B,...)alsoreturnsboththeouterand
%inneriterationnumbersatwhichXwascomputed:
0<
=ITER
(1)<
=MAXIT
%and0<
=ITER
(2)<
=RESTART.
%[X,FLAG,RELRES,ITER,RESVEC]=GMRES(A,B,...)alsoreturnsavectorof
%theresidualnormsateachinneriteration,includingNORM(B-A*X0).
%NotewithpreconditionersM1,M2,theresidualisNORM(M2\(M1\(B-A*X))).
%Example:
%n=21;
A=gallery('
wilk'
n);
b=sum(A,2);
%tol=1e-12;
maxit=15;
M=diag([10:
-1:
111:
10]);
%x=gmres(A,b,10,tol,maxit,M);
%Or,usethismatrix-vectorproductfunction
%%-----------------------------------------------------------------%
%functiony=afun(x,n)
%y=[0;
x(1:
n-1)]+[((n-1)/2:
0)'
;
(1:
(n-1)/2)'
].*x+[x(2:
n);
0];
%andthispreconditionerbacksolvefunction
%%------------------------------------------%
%functiony=mfun(r,n)
%y=r./[((n-1)/2:
1)'
1;
];
%asinputstoGMRES:
%x1=gmres(@(x)afun(x,n),b,10,tol,maxit,@(x)mfun(x,n));
%ClasssupportforinputsA,B,M1,M2,X0andtheoutputofAFUN:
%float:
double
%SeealsoBICG,BICGSTAB,BICGSTABL,CGS,LSQR,MINRES,PCG,QMR,SYMMLQ,
%TFQMR,ILU,FUNCTION_HANDLE.
%References
%H.F.Walker,"
ImplementationoftheGMRESMethodUsingHouseholder
%Transformations"
SIAMJ.Sci.Comp.Vol9.No1.January1988.
%Copyright1984-2010TheMathWorks,Inc.
%$Revision:
1.21.4.13$$Date:
2010/02/2508:
11:
48$
if(nargin<
2)
error('
MATLAB:
gmres:
NumInputs'
'
Notenoughinputarguments.'
);
end
%DeterminewhetherAisamatrixorafunction.
[atype,afun,afcnstr]=iterchk(A);
ifstrcmp(atype,'
matrix'
)
%Checkmatrixandrighthandsidevectorinputshaveappropriatesizes
[m,n]=size(A);
if(m~=n)
SquareMatrix'
Matrixmustbesquare.'
end
if~isequal(size(b),[m,1])
VectorSize'
'
%s%d%s'
...
'
Righthandsidemustbeacolumnvectoroflength'
m,...
tomatchthecoefficientmatrix.'
else
m=size(b,1);
n=m;
if~iscolumn(b)
Vector'
Righthandsidemustbeacolumnvector.'
%Assigndefaultvaluestounspecifiedparameters
3)||isempty(restart)||(restart==n)
restarted=false;
restarted=true;
4)||isempty(tol)
tol=1e-6;
warned=0;
iftol<
eps
warning('
tooSmallTolerance'
strcat('
Inputtolissmallerthanepsandmaynotbeachieved'
...
byGMRES\n'
Trytouseabiggertolerance'
));
warned=1;
tol=eps;
elseiftol>
=1
tooBigTolerance'
Inputtolisbiggerthan1\n'
Trytouseasmallertolerance'
tol=1-eps;
5)||isempty(maxit)
ifrestarted
maxit=min(ceil(n/restart),10);
else
maxit=min(n,10);
ifrestarted
outer=maxit;
ifrestart>
n
tooManyInnerIts'
RESTARTis%d%s\n%s%d.'
restart,'
butitshouldbeboundedbySIZE(A,1).'
SettingRESTARTto'
n);
restart=n;
inner=restart;
outer=1;
ifmaxit>
MAXITis%d%s\n%s%d.'
maxit,'
SettingMAXITto'
maxit=n;
inner=maxit;
%Checkforallzerorighthandsidevector=>
allzerosolution
n2b=norm(b);
%Normofrhsvector,b
if(n2b==0)%ifrhsvectorisallzeros
x=zeros(n,1);
%thensolutionisallzeros
flag=0;
%avalidsolutionhasbeenobtained
relres=0;
%therelativeresidualisactually0/0
iter=[00];
%noiterationsneedbeperformed
resvec=0;
%resvec
(1)=norm(b-A*x)=norm(0)
if(nargout<
itermsg('
gmres'
tol,maxit,0,flag,iter,NaN);
return
if((nargin>
=6)&
&
~isempty(M1))
existM1=1;
[m1type,m1fun,m1fcnstr]=iterchk(M1);
ifstrcmp(m1type,'
if~isequal(size(M1),[m,m])
PreConditioner1Size'
Preconditionermustbeasquarematrixofsize'
m,...
tomatchtheproblemsize.'
existM1=0;
m1type='
=7)&
~isempty(M2))
existM2=1;
[m2type,m2fun,m2fcnstr]=iterchk(M2);
ifstrcmp(m2type,'
if~isequal(size(M2),[m,m])
PreConditioner2Size'
existM2=0;
m2type='
=8)&
~isempty(x))
if~isequal(size(x),[n,1])
XoSize'
Initialguessmustbeacolumnvectoroflength'
n,...
8)&
strcmp(atype,'
)&
...
strcmp(m1type,'
strcmp(m2type,'
))
TooManyInputs'
Toomanyinputarguments.'
%Setupforthemethod
flag=1;
xmin=x;
%Iteratewhichhasminimalresidualsofar
imin=0;
%"
Outer"
iterationatwhichxminwascomputed
jmin=0;
Inner"
tolb=tol*n2b;
%Relativetolerance
evalxm=0;
stag=0;
moresteps=0;
maxmsteps=min([floor(n/50),5,n-maxit]);
maxstagsteps=3;
minupdated=0;
x0iszero=(norm(x)==0);
r=b-iterapp('
mtimes'
afun,atype,afcnstr,x,varargin{:
});
normr=norm(r);
%Normofinitialresidual
if(normr<
=tolb)%Initialguessisagoodenoughsolution
relres=normr/n2b;
resvec=normr;
tol,maxit,[00],flag,iter,relres);
minv_b=b;
ifexistM1
r=iterapp('
mldivide'
m1fun,m1type,m1fcnstr,r,varargin{:
if~x0iszero
minv_b=iterapp('
m1fun,m1type,m1fcnstr,b,varargin{:
minv_b=r;
if~all(isfinite(r))||~all(isfinite(minv_b))
flag=2;
x=xmin;
ifexistM2
m2fun,m2type,m2fcnstr,r,varargin{:
m2fun,m2type,m2fcnstr,minv_b,varargin{:
%normofthepreconditionedresidual
n2minv_b=norm(minv_b);
%normofthepreconditionedrhs
clearminv_b;
tolb=tol*n2minv_b;
relres=normr/n2minv_b;
resvec=n2minv_b;
resvec=zeros(inner*outer+1,1);
%Preallocatevectorfornormofresiduals
resvec
(1)=normr;
%resvec
(1)=norm(b-A*x0)
normrmin=normr;
%Normofresidualfromxmin
%PreallocateJtoholdtheGiven'
srotationconstants.
J=zeros(2,inner);
U=zeros(n,inner);
R=zeros(inner,inner);
w=zeros(inner+1,1);
foroutiter=1:
outer
%ConstructuforHouseholderreflector.
%u=r+sign(r
(1))*||r||*e1
u=r;
normr=norm(r);
beta=scalarsign(r
(1))*normr;
u
(1)=u
(1)+beta;
u=u/norm(u);
U(:
1)=u;
%ApplyHouseholderprojectiontor.
%w=r-2*u*u'
*r;
w
(1)=-beta;
foriniter=1:
inner
%FormP1*P2*P3...Pj*ej.
%v=Pj*ej=ej-2*u*u'
*ej
v=-2*(u(initer)'
)*u;
v(initer)=v(initer)+1;
%v=P1*P2*...Pjm1*(Pj*ej)
fork=(initer-1):
1
v=v-U(:
k)*(2*(U(:
k)'
*v));
%Explicitlynormalizevtoreducetheeffectsofround-off.
v=v/norm(v);
%ApplyAt
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- gmres 源程序 matlab