数值分析各种代码

追赶法functionx=tridiagsolver(a,b)[n,n]=size(a);fori=1:nifi==1l(i)=a(i,i);y(i)=b(i)/l(i);elsel(i)=a(i,i)-a(i,i-1)*u(i-1);y(i)=(b(i)-y(i-1)*a(i,i-1))/l(i);endifinu(i)=a(i,i+1)/l(i);endendx(n)=y(n);forj=n-1:-1:1x(j)=y(j)-u(j)*x(j+1);end拉格朗日functionyh=lagrange(x,y,xh)n=length(x);m=length(xh);yh=zeros(1,m);c1=ones(n-1,1);c2=ones(1,m);fori=1:nxp=x([1:i-1i+1:n]);yh=yh+y(i)*prod((c1*xh-xp'*c2)./(x(i)-xp'*c2));end线性x=[x1,x2]y=[y1.y2]xh=[xh]抛物线X=[x1,x2,x3]y=[y1，y2，y3]xh=[xh]牛顿差商（输入x,y为列向量）function[p,q]=chashang(x,y)n=length(x);p(:,1)=x;p(:,2)=y;forj=3:n+1p(1:n+2-j,j)=diff(p(1:n+3-j,j-1))./(x(j-1:n)-x(1:n+2-j));endq=p(1,2:n+1)';三次样条x=[0123];y=[0.200.52.01.5-1];pp=csape(x,y,'complete')[breaks,coefs,npolys,ncoefs,dim]=unmkpp(pp)最小二乘x=[0.240.650.951.241.732.012.232.522.772.99]';y=[0.23-0.26-1.1-0.450.270.1-0.290.240.561]';A=[log(x)cos(x)exp(x)];Z=A\y;a0=Z(1)a1=Z(2)a2=Z(3)x=[00.250.500.751.00];y=[1.001.2841.64872.11702.7183];p=polyfit(x,y,2)a2=p(1)a1=p(2)a0=p(3)复合中点functionI=fmid(fun,a,b,n)h=(b-a)/n;x=linspace(a+h/2,b-h/2,n);y=feval(fun,x);I=h*sum(y);复合梯形functionI=ftrapz(fun,a,b,n)h=(b-a)/n;x=linspace(a,b,n+1);y=feval(fun,x);I=h*(0.5*y(1)+0.5*y(n+1)+sum(y(2:n)));复合辛普森functionI=fsimpson(fun,a,b,n)h=(b-a)/n;x=linspace(a,b,2*n+1);y=feval(fun,x);I=h/6*(y(1)+y(2*n+1)+2*sum(y(3:2:2*n-1))+4*sum(y(2:2:2*n)));雅克比迭代function[x,iter]=jacobi(A,b,tol)D=diag(diag(A));L=D-tril(A);U=D-triu(A);x=zeros(size(b));foriter=1:500x=D\(b+U*x+L*x);error=norm(b-A*x)/norm(b);if(errortol)break;endendGS迭代function[x,iter]=GS(A,b,tol)D=diag(diag(A));L=D-tril(A);U=D-triu(A);x=zeros(size(b));foriter=1:500x=(D-L)\(b+U*x);error=norm(b-A*x)/norm(b);if(errortol)break;endendSOR迭代function[x,iter]=SOR(A,b,omega,tol)D=diag(diag(A));L=D-tril(A);U=D-triu(A);x=zeros(size(b));foriter=1:500x=(D-omega*L)\(omega*b+(1-omega)*D*x+omega*U*x);error=norm(b-A*x)/norm(b);if(errortol)break;endend二分法functionv=f(x)v=x^3-x-1;functionx=erfenfa(f,a,b,tol)ifnargin4tol=1e-5;endfa=feval(f,a);fb=feval(f,b);whileabs(a-b)tolx=(a+b)/2;fx=feval(f,x);ifsign(fx)==sign(fa)a=x;fa=fx;elseifsign(fx)==sign(fb)b=x;fb=fx;endendendx=erfenfa('f',1,2)牛顿法function[x,it]=newton(f,g,x0,tol)it=0;done=0;while~donex=x0-feval(f,x0)/feval(g,x0);it=it+1;done=(norm(x-x0)=tol);if~donex0=x;endendfunctionr=f(x)r=polyval([1,2,10,-20],x);functionr=g(x)p=polyder([1,2,10,-20]);r=polyval(p,x);牛顿下山法乘幂法function[t,y]=chengmifa(a,xinit,ep)v0=xinit;[tv,ti]=max(abs(v0));lam0=v0(ti);u0=v0/lam0;flag=0;while~flagv1=a*u0;[tv,ti]=max(abs(v1));lam1=v1(ti);u0=v1/lam1;err=abs(lam0-lam1);iferrepflag=1endlam0=lam1;endt=lam1;y=u0;反幂法function[t,y]=fanmifa(a,xinit,ep)v0=xinit;[tv,ti]=max(abs(v0));lam0=v0(ti);u0=v0/lam0;flag=0;while~flagv1=inv(a)*u0;[tv,ti]=max(abs(v1));lam1=v1(ti);u0=v1/lam1;err=abs(lam0-lam1);iferrepflag=1;endlam0=lam1;endt=1/lam1;y=u0;改进欧拉法function[x,y]=odeIEuler(f,y0,a,b,n)h=(b-a)/n;x=a:h:b;y(1)=y0;fori=1:nyp=y(i)+h*feval(f,x(i),y(i));yc=y(i)+h*feval(f,x(i+1),yp);y(i+1)=1/2*(yp+yc);endfunctionr=f(x,y)r=y-2*x/y;