数值分析的几个简单算法实现.docVIP

数值分析的几个简单算法实现??Matlab6.5 ? ? 2006-11-24 07:42:12 大 中 小 %埃特金算法Aitken.m function x=Aitken(fname,x0,N,e) k=1;e=1e-5;x2=x0+2*e; while abs(x2-x0)ekN x1=feval(x0);x2=feval(x1); x2=x2-(x2-x1)^2/(x2-2*x1+x0); k=k+1; x2=x0;end; disp(x2); if k==N,warning(迭代失败);end; end ? ? ? %迭代法diedai.m文件 function x=diedai(fname,x0,e,N) if nargin5;N500;end; if nargin4;e=1e-4;end; x=x0;x0=x+2*e;e=1e-5;k=1; while abs(x0-x)ekN; x=feval(fname,x0); k=k+1; x0=x; disp(x); end if k==N,warning(迭代失败);end ? ? %二分法erfenfa.m function x=erfenfa(fname,a,b,e) y0=feval(fname,a);e=1e-5; while (b-a)e x=(a+b)/2; y=feval(fname,x); if y*y00 a=x; else b=x; end disp(x); end ? ?%牛顿法 function x=newton(fname,dfname,x0,e,N) if nargin5;N500;end; if nargin4;e=1e-4;end; x=x0;x0=x+2*e;k=0; while abs(x0-x)ekN k=k+1; x0=x;x=x0-feval(fname,x0)/feval(dfname,x0); disp(x); end if k==N,warning(已达迭代次数上限);end % format long; newton(inline(exp(x)-4*cos(x)),inline(exp(x)+4*sin(x)),pi/4,1e-4,20) ? ? ? ? function t=Rombeg(fname,a,b,e) %用途：龙贝格法求函数积分 %格式：t=Rombeg(fname,a,b,e) fname是被积函数，a,b分别为下上限，e为精度（默认1e-4） if nargin4,e=1e-4;end; i=1;j=1;h=b-a; T(i,1)=h/2*(feval(fname,a)+feval(fname,b)); T(i+1,1)=T(i,1)/2+sum(feval(fname,a+h/2:h:1-h/2+0.001*h))*h/2; T(i+1,j+1)=4^j*T(i+1,j)/(4^j-1)-T(i,j)/(4^j-1); while abs(T(i+1,i+1)-T(i,i))e i=i+1;h=h/2; T(i+1,1)=T(i,1)/2+sum(feval(fname,a+h/2:h:1-h/2+0.001*h))*h/2; for j=1:i T(i+1,j+1)=4^j*T(i+1,j)/(4^j-1)-T(i,j)/(4^j-1); end end T t=T(i+1,j+1); ? % format long;Rombeg(inline(sin(x)./x),eps,1,1e-6) function x=gaussseidel(A,b,x0,e,N) % 用途：用向量形式的gauss-seidel迭代解线性方程组Ax=b if nargin5,N=500;end x=x0;x0=x+2*e; k=0;A1=tril(A);A2=inv(A1); while norm(x0-x,inf)ekN, k=k+1; x0=x;x=-A2*(A-A1)*x0+A2*b; disp(x) end if k==N,warning(已达到迭代次数上限); end x=gaussseidel(A,b,[1 1 1],1e-6) ? ? function x=sor(A,b,omega,x0,e,N) % 用途：用分量形式的sor迭代解线性方程组Ax=b n=length(b); if nargin6,N=500;end x=x0;x0=x+2*e; k=0;L=