matlab实现插值法和曲线拟合

1、插值法和曲线拟合电子科技大学 摘要：理解拉格朗日多项式插值、分段线性插值、牛顿前插，曲线拟合，用matlab编程求解函数，用插值法和分段线性插值求解同一函数，比较插值余项；用牛顿前插公式计算函数，计算函数值；对于曲线拟合，用不同曲线拟合数据。关键字：拉格朗日插值多项式；分段线性插值；牛顿前插；曲线拟合引言：在数学物理方程中，当给定数据是不同散点时，无法确定函数表达式，求解函数就需要很大的计算量，我们有多种方法对给定的表格函数进行求解，我们这里，利用插值法和曲线拟合对函数进行求解，进一步了解函数性质，两种方法各有利弊，适合我们进行不同的散点函数求解。正文：一、插值法和分段线性插值1拉格朗日多项式

2、原理对某个多项式函数，已知有给定的k + 1个取值点：其中对应着自变量的位置，而对应着函数在这个位置的取值。假设任意两个不同的xj都互不相同，那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为：其中每个为拉格朗日基本多项式（或称插值基函数），其表达式为：3拉格朗日基本多项式的特点是在 上取值为1，在其它的点 上取值为0。2分段线性插值原理给定区间a,b, 将其分割成a=x0 x1 =x0(i) & xi(j)=x0(i+1) & k=3 lx(1)=(xi(j)-x0(i+1)/(x0(i)-x0(i+1); lx(2)=(xi(j)-x0(i)/(x0(i+1)-x0(i); y(k)=l

3、x(1)*y0(i)+lx(2)*y0(i+1); k=k+1; end endendplot(xi,y,r);axis(0.2 0.8 1.03 1.24);hold on plot(x0,y0,b.,markersize,20)grid on6运算结果拉格朗日插值结果x=0.,y=1.x=0.,y=1.x=0.,y=1.拉格朗日插值余项：分段插值结果ans =1.0508 1.1418 1.2095分段线性插值余项：由于拉格朗日插值的余项比分段线性插值的余项要求更为严格，点少、区间小的时候，拉格朗日插值要更好。但在区间较大、节点较多的时候，分段线性插值要更好。二、牛顿前插1牛顿前插原理次牛

4、顿前插公式：插值余项：，阶差分记作。 阶差商是差分和差商之间的关系是2牛顿前插算法输入。对，计算各阶差分计算函数值3牛顿前插程序：编写一个用牛顿前插公式计算函数值的程序，要求先输出差分表，再计算x点的函数值0.1250.2500.3750.5000.6250.7500.7960.7730.7440.7040.6560.602分别求x=0.158和x=0.636的三次插值的值，并比较二者的插值余项。这里以x=0.636为例function P=newtonchax0=0.636;X=0.125 0.250 0.375 0.500 0.625 0.750;Y=0.796 0.773 0.744 0

5、.704 0.656 0.602;h=abs(X(2)-X(1);n=find(abs(x0-X)3*h);X=X(n(1):n(end);Y=Y(n(1):n(end);w=length(X);R=zeros(w,w);R(:,1)=Y(:);for k=2:w for j=k:w R(j,k)=R(j,k-1)-R(j-1,k-1); endendt=(x0-X(1)/h;T=1; for m=1:w-1 T=T*(t-m+1); N(m)=R(m+1,m+1)*T/factorial(m); endP=R(1,1)+sum(N);4运行结果：差分表0.0000000000.0000-0.

6、000000000.0000-0.0000-0.00000000.0000-0.0000-0.0000-0.0000000.0000-0.0000-0.00000.00000.000000.0000-0.0000-0.00000.0000-0.0000-0.0000X=0.636时ans =0.4000 x=0,158时ans =0.0000三、曲线拟合1曲线拟合原理：给定数据。记拟合函数的形式为（1.1），其中为已知的线性无关函数。求系数使得（1.2）取最小值。称（1.3）为拟合函数或经验公式。如果，则（1.3）为次最小二乘拟合多项式2曲线拟合算法：已知数据对，求多项式，使得为最小。注意到此

7、时，多项式系数满足下面的线性方程组：其中，然后只要调用线性方程组的函数程序即可3曲线拟合程序：试分别用抛物线y=a+bx2和指数曲线y=aebx拟合下列数据12.53.543.81.5026.033.0画出数据点和两条拟合曲线，并通过计算2个拟合函数残差向量的2范数来比较拟合优劣。用抛物线y=a+bx拟合程序：function ZXEx=1 2.52 3.52 42;y=3.8 1.50 26.0 33.0;m=1;S=zeros(1,2*m+1);T=zeros(m+1,1);for k=1:2*m+1 S(k)=sum(x.(k-1);endfor k=1:m+1 T(k)=sum(x.(

8、k-1).*y);endA=zeros(m+1,m+1);a=zeros(m+1,1);for i=1:m+1 for j=1:m+1 A(i,j)=S(i+j-1); endenda=AT;for k=1:m+1 fprintf(a%d=%fn,k,a(k);endp=polyfit(x,y,1);u=polyval(p,x);plot(sqrt(x),u,b)hold onplot(sqrt(x),y,b.)grid on指数曲线y=aebx拟合程序：function ZXE2x=1 2.5 3.5 4;y=3.8 1.50 26.0 33.0;y=log(y);m=1;S=zeros(1,2*m+1);T=zeros(m+1,1);for k=1:2*m+1 S(k)=sum(x.(k-1);endfor k=1:m+1 T(k)=sum(x.(k-1).*y);endA=zeros(m+1,m+1);a=zeros(m+1,1);for i=1:m+1 for j=1:m+1 A(i,j)=S(i+j-1); endenda=AT;for k=1:m+1 fprintf(a%d=%fn,k,a(k);endp=polyfit(x,y,1);u=polyval(p,x);plot(x,exp(u),r)hold onplot(x,ex