雅可比迭代法:function x=jacobi(a,b,p,delta,n)%a为n维非奇异矩阵;b为n维值向量%p为初值;delta为误差界;n为给定的迭代最高次数N=length(b);for k=1:nfor j=1:Nx(j)=(b(j)-a(j,[1:j-1,j+1:N])*p([1:j-1,j+1:N]))/a(j,j);enderr=abs(norm(x’-p));p=x’;if(err<delta)break;endendp %显示迭代过程x=x’;k,err高斯塞德尔法迭代:function x=saidel(a,b,p,delta,n)%a为n维非奇异矩阵;b为n维值向量%p为初值;delta为误差界;n为给定的迭代最高次数N=length(b);for k=1:nfor j=1:Nif j==1x(1)=(b(1)-a(1,2:N)*p(2:N))/a(1,1);else if j=Nx(N)=(b(N)-a(N,1:N-1)*(x(1:N-1))’)/a(N,N);elsex(j)=(b(j)-a(j,(1:j-1)*x(1:j-1)-a(j,j+1:N)*p(j+1:N))/a(j,j);endenderr=abs(norm(x’-p));p=x’;if(err<delta)break;endendx=x’;k,err不动点迭代法:function [x,k,err,p]=ddf(f,x0,tol,n)%ddl.m为用迭代法求非线性方程的解%f为给定的迭代函数;x0为给定的初始值%tol为给定的误差界;n为所允许的最大迭代次数%k为迭代次数;x为不动点的近似值;err为误差p(1)=x0;for k=2:np(k)=feval(f,p(k-1));k,err=abs(p(k)-p(k-1))x=p(k);if(err<tol)break;endif k==ndisp('迭代超过最大次数!')endendx=p'牛顿法:function [x,k,err,y]=Newtun(f,df,x0,tol,n)%Newtun.m为用迭代法求非线性方程的解%f为给定的非线性方程;df为f的微分方程;x0为给定的初始值%tol为给定的误差界;n为所允许的最大迭代次数%k为迭代次数;x为不动点的近似值;err为误差%x为牛顿迭代法得到得近似解y(1)=x0;for k=1:nx=x0-feval(f,x0)/feval(df,x0);err=abs(x-x0);x0=x;if(err<tol)|(y==0);break;endend必要编辑M文件qfun.m,代码如下:function y=qfun(x);y=x^3-3*x-1;弦截法:function [x,err,k,y]=xjf(f,x0,x1,tol,n)%xjf.m为用弦截法迭代法求非线性方程的解%f为给定的非线性方程;x0,x1为给定的初始值%tol为给定的误差界;n为所允许的最大迭代次数%k为迭代次数;x为牛顿迭代法的近似值;err为x1-x0的绝对值y(1)=x0;for k=1:nx=x1-feval('li6_5fun',x1)*(x1-x0)/(feval('li6_5fun',x1)-feval('li6_5fun',x0));err=abs(x-x1);x0=x1;x1=x;if(err<tol)|(y==0);break;endend必要编辑M文件li6_5.m,代码如下:function y=li6_5(x);y=x^3-3*x-1;复化梯形公式matlab:function t=tixing(f,a,b,n)h=(b-a)/n;sum=0;for k=0:n-1x=a+k*h;sum=sum+feval(f,x);endt=h/2*(feval(f,a)+feval(f,b)+2*sum);运行程序结果:>> format long>> tixing(inline('x./(x.^2+4)'),0,1,8)ans =0.111402354529548复化辛普森公式matlab:function y=xinpusen(f,a,b,n)h=(b-a)/n;sn=h/6*(feval(f,a)+feval(f,b));s=0;for i=0:n-1s=s+4*feval(f,a+h/2+i*h)+2*feval(f,a+i*h);endy=sn+h*s/6;运行结果:>> xinpusen(inline('x./(x.^2+4)'),0,1,8)ans =0.111571813252631改进欧拉法:function y=gjEuler(f,x0,y0,xn,N)%gjEuler.m函数为改进的欧拉法求微分方程的解%f为一阶常微分方程的一般表达式的右端函数%x0,y0为初始条件%xn为取值范围的一个端点%N为区间的个数%y为求解微分方程组的解x=zeros(1,N+1);y=zeros(1,N+1);x(1)=x0;y(1)=y0;h=(xn-x0)/N;for n=1:Nx(n+1)=x(n)+h;z0=y(n)+h*feval(f,x(n),y(n));y(n+1)=y(n)+h/2*(feval(f,x(n),y(n))+feval(f,x(n+1),z0));endT=[x',y']梯形法:function y=TXF(f,x0,y0,xn,N)%TXF.m函数为改进的梯形法求微分方程的解%f为一阶常微分方程的一般表达式的右端函数%x0,y0为初始条件%xn为取值范围的一个端点%N为区间的个数%y为求解微分方程组的解x=zeros(1,N+1);y=zeros(1,N+1);x(1)=x0;y(1)=y0;h=(xn-x0)/N;for n=1:Nx(n+1)=x(n)+h;y(n+1)=y(n)+h/2*(feval(f,x(n),y(n))+feval(f,x(n+1),y(n+1)));endT=[x',y']function y=LG4(f,x0,y0,xn,N)%LG4.m函数为四阶龙格库塔法求微分方程的解%f为一阶微分方程%x0,y0为左右端点%xn为给定的初始值%N为给定的迭代步长%y为求解微分方程组的解x=zeros(1,N+1);y=zeros(1,N+1);h=(y0-x0)/N;T=x0:h:y0;y(1)=xn;for n=1:Nk1=h*feval(f,T(n),y(n));k2=h*feval(f,T(n)+h/2,y(n)+k1/2);k3=h*feval(f,T(n)+h/2,y(n)+k2/2);k4=h*feval(f,T(n)+h,y(n)+k3);y(n+1)=y(n)+(k1+2*k2+2*k3+k4)/6;endR=[T',y']直接三角分解:function dxl=dxllu(a,b)format rat[n n] =size(a);RA=rank(a);if RA~=ndisp('请注意:因为A的n阶行列式hl等于零,所以A不能进行LU分解.A的秩RA如下:')RAhl=det(a)returnendif RA==nfor p=1:nh(p)=det(a(1:p, 1:p));endhl=h(1:n);for i=1:nif h(i)==0disp('请注意:因为A的r阶主子式等于零,所以A不能进行LU分解.A的秩RA 和各阶顺序主子式值hl依次如下:')hlRAreturnendenddisp('请注意:因为A的各阶主子式都不等于零,所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下:')u(1,1:n)=a(1,1:n);l(2:n,1)=a(2:n,1)/u(1,1);for i=1:nl(i,i)=1;endfor k=2:nfor j=k:ns=0;for m=1:(k-1)s=s+l(k,m)*u(m,j);endu(k,j)=a(k,j)-s;endfor i=(k+1):ns=0;for m=1:(k-1)s=s+l(i,m)*u(m,k);endl(i,k)=(a(i,k)-s)/u(k,k);endendy(1)=b(1);for i=2:np=0;for k=1:i-1p=p+l(i,k)*y(k);endy(i)=b(i)-p;endx(n)=y(n)/u(n,n);for i=(n-1):(-1):1q=0;for k=i+1:nq=q+u(i,k)*x(k);endx(i)=(y(i)-q)/u(i,i);enddisp('各阶主子式:')hldisp('系数矩阵的秩:')RAdisp('单位下三角:') ldisp('上三角:')udisp('结果Y为:') ydisp('结果x为:') xend。