第一章插值方法1.1Lagrange插值1.2逐步插值1.3分段三次Hermite插值1.4分段三次样条插值第二章数值积分2.1 Simpson公式2.2 变步长梯形法2.3 Romberg加速算法2.4 三点Gauss公式第三章常微分方程德差分方法3.1 改进的Euler方法3.2 四阶Runge-Kutta方法3.3 二阶Adams预报校正系统3.4 改进的四阶Adams预报校正系统第四章方程求根4.1 二分法4.2 开方法4.3 Newton下山法4.4 快速弦截法第五章线性方程组的迭代法5.1 Jacobi迭代5.2 Gauss-Seidel迭代5.3 超松弛迭代5.4 对称超松弛迭代第六章线性方程组的直接法6.1 追赶法6.2 Cholesky方法6.3 矩阵分解方法6.4 Gauss列主元消去法第一章插值方法1.1Lagrange插值计算Lagrange插值多项式在x=x0处的值. MATLAB文件：（文件名：Lagrange_eval.m）function [y0,N]= Lagrange_eval(X,Y,x0)%X,Y是已知插值点坐标%x0是插值点%y0是Lagrange插值多项式在x0处的值%N是Lagrange插值函数的权系数m=length(X);N=zeros(m,1);y0=0;for i=1:mN(i)=1;for j=1:mif j~=i;N(i)=N(i)*(x0-X(j))/(X(i)-X(j));endendy0=y0+Y(i)*N(i);end用法》X=[…];Y=[…];》x0= ;》[y0,N]= Lagrange_eval(X,Y,x0)1.2逐步插值计算逐步插值多项式在x=x0处的值.MATLAB文件：（文件名：Neville_eval.m）function y0=Neville_eval(X,Y,x0)%X,Y是已知插值点坐标%x0是插值点%y0是Neville逐步插值多项式在x0处的值m=length(X);P=zeros(m,1);P1=zeros(m,1);P=Y;for i=1:mP1=P;k=1;for j=i+1:mk=k+1;P(j)=P1(j-1)+(P1(j)-P1(j-1))*(x0-X(k-1))/(X(j)-X(k-1));endif abs(P(m)-P(m-1))<10^-6;y0=P(m);return;endendy0=P(m);用法》X=[…];Y=[…];》x0= ;》y0= Neville_eval(X,Y,x0)1.3 分段三次Hermite插值利用分段三次Hermite插值计算插值点处的函数近似值. MATLAB文件：（文件名：Hermite_interp.m）function y0=Hermite_interp(X,Y,DY,x0)%X,Y是已知插值点向量序列%DY是插值点处的导数值%x0插值点横坐标%y0为待求的分段三次Hermite插值多项式在x0处的值%N向量长度N=length(X);for i=1:Nif x0>=X(i)&& x0<=X(i+1)k=i; break;endenda1=x0-X(k+1);a2=x0-X(k);a3= X(k)- X(k+1);y0=(a1/a3)^2*(1-2*a2/a3)*Y(k)+(-a2/a3)^2*(1+2*a1/a3)*Y(k+1)+...(a1/a3)^2*a2*DY(k)+(-a2/a3)^2*a1*DY(k+1);用法》X=[…];Y=关于X的函数;DY=Y’;》x0=插值横坐标;》y0==Hermite_interp(X,Y,DY,x0)1.4分段三次样条插值计算在插值点处的函数值，并用来拟合曲线.MATLAB文件：（文件名：Spline_interp.m）function [y0,C]=Spline_interp(X,Y,s0,sN,x0)%X,Y是已知插值坐标%s0,sN是两端点的一次导数值%x0是插值点%y0三次样条函数在x0处的值%C是分段三次样条函数的系数N=length(X);C=zeros(4,N-1); h=zeros(1,N-1);mu=zeros(1,N-1); lmt=zeros(1,N-1);d=zeros(1,N); %d表示右端函数值h=X(1,2:N)-X(1,1:N-1);mu(1,N-1)=1; lmt(1,1)=1;mu(1,1:N-2)=h(1,1:N-2)/(h(1,1:N-2)+h(1,2:N-1));lmt(1,2:N-1)=h(1,2:N-1)/(h(1,1:N-2)+h(1,2:N-1));d(1,1)=6*((Y(1,2)-Y(1,1))/h(1,1)-s0)/h(1,1);d(1,N)=6*(sN-(Y(1,N)-Y(1,N-1))/h(1,N-1))/h(1,N-1);d(1,2:N-1)=6*((Y(1,3:N)-Y(1,2:N-1))/h(1,2:N-1)…(Y(1,2:N-1)-Y(1,1:N-2))/h(1,1:N-2))/(h(1,1:N-2)+h(1,2:N-1)); %追赶法解三对角方程组bit=zeros(1,N-1);bit(1,1)=lmt(1,1)/2;for i=2:N-1bit(1,i)=lmt(1,i)/(2-mu(1,i-1)*bit(1,i-1));endy=zeros(1,N);y(1,1)=d(1,1)/2;for i=2:Ny(1,i)=(d(1,i)-mu(1,i-1)*y(1,i-1))/(2-mu(1,i-1)*bit(1,i-1)); endx=zeros(1,N);x(1,N)=y(1,N);for i=N-1:-1:1x(1,i)=y(1,i)-bit(1,i)*x(1,i+1);endv=zeros(1,N-1);v(1,1:N-1)=(Y(1,2:N)-Y(1,1:N-1))/h(1,1:N-1);C(4,:)=Y(1,1:N-1);C(3,:)=v-h*(2*x(1,1:N-1)+x(1,2:N))/6;C(2,:)=x(1,1:N-1)/2;C(1,:)=(x(1,2:N)-x(1,1:N-1))/(6*h);if nargin<5y0=0;elsefor j=1:N-1if x0>=X(1,j)& x0<X(1,j+1)omg=x0-X(1,j);y0=((C(4,j)*omg+C(3,j))omg+C(2,j))*omg+C(1,j);endendend试求三次样条插值函数S(x)，其中：S’(0.25)=1.0000,S’(0.53)=0.6868. 解：X=[0.25,0.30,0.39,0.45,0.53];Y=[0.5000,0.5477,0.6254,0.6708,0.7280]; s0=1.0000;sN=0.6868;[y0,C]=Spline_interp(X,Y,s0,sN,x0);plot(0.25:0.01:0.30,polyval(C(:.1),0:0.01:0.05),’r-.’);hold onplot(0.30:0.01:0.39,polyval(C(:.2),0:0.01:0.09),’b’);plot(0.39:0.01:0.45,polyval(C(:.3),0:0.01:0.06),’k-*’);plot(0.45:0.01:0.53,polyval(C(:.4),0:0.01:0.08));第二章数值积分2.1 Simpson公式利用复化Simpson公式求被积函数f(x)在给定区间上的积分值. MATLAB文件：（文件名：FSimpson.m）function S=FSimpson(f,a,b,N)%f表示被积函数句柄%a,b表示被积区间端点%N表示区间个数%S是用复化Simpson公式求得的积分值h=(b-a)/N;fa=feval(f,a);fb=feval(f,b);S=fa+fb;x=a;for i=1:Nx=x+h/2;fx=feval(f,x);S=S+4*fx;x=x+h/2;fx=feval(f,x);S=S+2*fx;endS=h*S/6;算例2：利用复化Simpson 公式计算积分⎰+=1024x xdx S . 解：后面都要用到的f1:function f=f1(x)f=x/(4+x^2);令f=@f1;a=0;b=1;运行S=FSimpson(f,a,b,N)这里的N 值需要自己输入。

2.2 变步长梯形法利用变步长梯形法求被积函数f(x)在给定区间上的积分值. MATLAB 文件：（文件名：bbct.m ）function [T,n]=bbct(f,a,b,eps)%f 表示被积函数句柄%a,b 被积区间端点%eps 精度%T 是用变步长梯形法求得的积分值%n 表示二分区间的次数h=b-a;fa=feval(f,a);fb=feval(f,b);T1=h*(fa+fb)/2;T2=T1/2+h*feval(f,a+h/2)/2;n=1;%按变步长梯形法求积分值while abs(T2-T1)>=epsh=h/2;T1=T2;S=0;x=a+h/2;while x<bfx=feval(f,x);S=S+fx;x=x+h;endT2=T1/2+S*h/2;n=n+1;endT=T2;算例3：利用变步长梯形法计算积分⎰+=1024x xdx T . 解：function f=f1(x)f=x/(4+x^2);令f=@f1;a=0;b=1;运行 [T,n]=bbct(f,a,b,eps)这里的eps 值需要自己输入。

2.3 Romberg 加速算法利用Romberg 加速算法计算被积函数f(x)在给定区间上的积分值. MATLAB 文件：（文件名：Romberg.m ）function [quad,R]=Romberg(f,a,b,eps)%f 表示被积函数句柄%a,b 被积区间端点%eps 精度%quad 是用Romberg 加速算法求得的积分值%R 为Romberg 表%err 表示误差的估计h=b-a;R(1,1)=h*(feval(f,a)+feval(f,b))/2;M=1; J=0; err=1;while err>epsJ=J+1;h=h/2;S=0;for p=1:Mx=a+h*(2*p-1);S=S+feval(f,x);endR(J+1,1)=R(J,1)/2+h*S;M=2*M;for k=1:JR(J+1,k+1)=R(J+1,k)+(R(J+1,k)-R(J,k))/(4^k-1);enderr=abs(R(J+1,J)-R(J+1,J+1));endquad=R(J+1,J+1);算例4：利用Romberg 加速算法计算积分⎰+=1024x xdx R . 解：function f=f1(x)f=x/(4+x^2);令f=@f1;a=0;b=1;运行 [quad,R]=Romberg(f,a,b,eps)这里的eps 值需要自己输入。