PAGEPAGE14PAGE1實(shí)驗(yàn)二插值與擬合實(shí)驗(yàn)名稱：插值與擬合實(shí)驗(yàn)類型:驗(yàn)證性實(shí)驗(yàn)學(xué)時(shí)：23.1實(shí)驗(yàn)環(huán)境操作系統(tǒng)：WindowsXP/Win7編程環(huán)境：自定3.2實(shí)驗(yàn)?zāi)康恼莆斩囗?xiàng)式插值法的基本思路和步驟；了解整體插值的局限性及分段插值的基本思想。掌握最小二乘法擬合的基本原理和方法；培養(yǎng)運(yùn)用計(jì)算機(jī)模擬解決問題的能力。3.3實(shí)驗(yàn)原理和方法3.3.1多項(xiàng)式插值若n次多項(xiàng)式在n+1個(gè)插值節(jié)點(diǎn)上滿足插值條件：則稱這n＋1個(gè)n次多項(xiàng)式為插值節(jié)點(diǎn)上的n次插值基函數(shù)。由于時(shí)，，故為的零點(diǎn)，從而可以設(shè)由可得所以可得。若記，則有，從而有。上式是次數(shù)不超過的多項(xiàng)式，且滿足所有的插值條件，稱之為L(zhǎng)agrange插值多項(xiàng)式。3.3.2數(shù)據(jù)擬合最小二乘法基本原理已知數(shù)據(jù)對(duì)，求多項(xiàng)式使得為最小，這就是多項(xiàng)式擬合的最小二乘法。 最小二乘法的算法描述 線性函數(shù)為例，擬合給定數(shù)據(jù)。算法描述：步驟1：輸入值，及。步驟2：建立法方程組，并求解。步驟3：輸出。3.4實(shí)驗(yàn)內(nèi)容和步驟3.4.1多項(xiàng)式插值1．給定構(gòu)造插值多項(xiàng)式計(jì)算。編程實(shí)現(xiàn)拉格朗日插值，并計(jì)算結(jié)果。將計(jì)算結(jié)果和查表結(jié)果進(jìn)行比較。區(qū)間作等距劃分：，以（）為節(jié)點(diǎn)對(duì)函數(shù)進(jìn)行插值逼近。（分別取）用多項(xiàng)式插值對(duì)進(jìn)行逼近，并在同一坐標(biāo)系下作出函數(shù)的圖形，進(jìn)行比較。寫出插值函數(shù)對(duì)的逼近程度與節(jié)點(diǎn)個(gè)數(shù)的關(guān)系，并分析原因。試用分段插值（任意選?。?duì)進(jìn)行逼近，在同一坐標(biāo)下畫出圖形，觀察分段插值函數(shù)對(duì)的逼近程度與節(jié)點(diǎn)個(gè)數(shù)的關(guān)系。

解：1.（1）拉格朗日插值程序代碼：functionf=Language(x,y,x0)symstl;if(length(x)==length(y))n=length(x);elsedisp('x和y的維數(shù)不相等！');return;%檢錯(cuò)endh=sym(0);for(i=1:n)l=sym(y(i));for(j=1:i-1)l=l*(t-x(j))/(x(i)-x(j));end;for(j=i+1:n)l=l*(t-x(j))/(x(i)-x(j));end;h=h+l;endsimplify(h);if(nargin==3)f=subs(h,'t',x0);%計(jì)算插值點(diǎn)的函數(shù)值elsef=collect(h);f=vpa(f,6);%將插值多項(xiàng)式的系數(shù)化成6位精度的小數(shù)endx=[11,12,13];y=[0.190809,0.207912,0.224951];f=Language(x,y,11.5);vpa(f,8)

程序結(jié)果：查表所得的精確解：0.199367934，所以插值所得到的結(jié)果精確到小數(shù)點(diǎn)后五位。

解：2.程序代碼：functiony=lagrange(x0,y0,x);n=length(x0);m=length(x);fori=1:mz=x(i);s=0.0;fork=1:np=1.0;forj=1:nifj~=kp=p*(z-x0(j))/(x0(k)-x0(j));endends=p*y0(k)+s;endy(i)=s;endx1=[-5:2:5];y=1./(1+x1.^2);x=[-5:0.001:5];y1=L(x1,y);x2=[-5:1:5];y=1./(1+x2.^2);y2=L(x2,y);x3=[-5:0.5:5];y=1./(1+x3.^2);y3=L(x3,y);y=1./(1+x3.^2);plot(x,subs(y1,'t',x),'r',x,subs(y2,'t',x),'g',x,subs(y3,'t',x),'b')holdonplot(x,1./(1+x.*x),'m')axis([-55-0.52])xlabel('x')ylabel('y')legend('n=5插值后的圖','n=10插值后的圖','n=20插值后的圖','原圖',0)結(jié)果與原因分析：通過上圖發(fā)現(xiàn)，高次插值穩(wěn)定性差，而低次插值對(duì)于較大區(qū)間逼近精度又不夠，而且，隨著插值節(jié)點(diǎn)數(shù)目的增多，區(qū)域中間的確越來(lái)越逼近原函數(shù)，然而區(qū)域兩端發(fā)生了較大震蕩，即發(fā)生了Runge現(xiàn)象，這是由叉子余項(xiàng)的表達(dá)式所決定的。解決這一矛盾的有效方法就是采用分段低次代數(shù)插值。分段插值函數(shù)程序代碼：functionoutput=L1symsfxlx;f=1/(1+x^2);N=input('請(qǐng)輸入插值節(jié)點(diǎn)數(shù)N=');xx=-5:10/N:5;ff=zeros(1,length(xx));fori=1:(N+1)x=xx(i);ff(i)=eval(f);endM=-5:0.01:5;output=zeros(1,length(M));n=1;fori=2:N+1forx=-5:0.01:5ifx<xx(i)&&x>=xx(i-1)lx(1)=ff(i-1)*(x-xx(i))/(xx(i-1)-xx(i));lx(2)=ff(i)*(x-xx(i-1))/(xx(i)-xx(i-1));output(n)=lx(1)+lx(2);n=n+1;endendendA=-5:0.01:5;Endoutput=L1;plot(A,output,'r');holdonoutput=L1;plot(A,output,'g');holdonoutput=L1;plot(A,output,'b');holdonezplot(f,[-5,5])xlabel('x')ylabel('y')legend('n=5分段差值后的圖','n=10分段差值后的圖','n=20分段插值后的圖','原圖',0)請(qǐng)輸入插值節(jié)點(diǎn)數(shù)N=5請(qǐng)輸入插值節(jié)點(diǎn)數(shù)N=10請(qǐng)輸入插值節(jié)點(diǎn)數(shù)N=20程序結(jié)果：結(jié)果分析：由圖可以發(fā)現(xiàn)，插值節(jié)點(diǎn)數(shù)目越多，插值函數(shù)越接近原函數(shù)，函數(shù)沒有出現(xiàn)Runge現(xiàn)象，這是因?yàn)榉侄尾逯祵^(qū)間分為了n份。3.4.2數(shù)據(jù)擬合1．已知一組數(shù)據(jù)如下，求它的線性擬合曲線。1234544.5688.5編程實(shí)現(xiàn)最小二乘算法，并畫出其擬合曲線求出其平方誤差。2．已知一組數(shù)據(jù)如下，求其擬合曲線。01234567891023478101114161819106.42108.2109.5110109.93110.49110.59110.6110.76111111.2求以上數(shù)據(jù)形如的擬合曲線，及其平方誤差。求以上數(shù)據(jù)形如的擬合曲線，及其平方誤差。通過畫出（1）（2）的圖形，觀察結(jié)果并結(jié)合其平方誤差，寫出你對(duì)數(shù)據(jù)擬合的認(rèn)識(shí)。解：1.最小二乘法程序代碼x=[12345];y=[4,4.5,6,8,8.5];x1=0;x12=0;y1=0;xy=0;n=5;fori=1:nx1=x1+x(i);x12=x12+x(i)^2;y1=y1+y(i);xy=xy+x(i)*y(i);enda0=(y1*x12-x1*xy)/(n*x12-x1*x1);a1=(n*xy-x1*y1)/(n*x12-x1*x1);yy=a0+a1*xplot(x,y,'*')holdonplot(yy,'r')axis([0,6,0,9])結(jié)果：（2）平方誤差：2.（1）MATLAB程序代碼：

xx=0;yy=0;x2=0;x3=0;x4=0;xy=0;y2=0;n=11;x1=[2,3,4,7,8,10,11,14,16,18,19];y1=[106.42,108.2,109.5,110,109.93,110.49,110.59,110.6,110.76,111,111.2];fori=1:nxx=xx+x1(i);x2=x2+x1(i)^2;x3=x3+x1(i)^3;x4=x4+x1(i)^4;yy=yy+y1(i);xy=xy+x1(i)*y1(i);y2=y2+x1(i)*x1(i)*y1(i);endA=[nxxx2;xxx2x3;x2x3x4];B=[yy;xy;y2];C=A\B;y2=C(1)+C(2)*x1+C(3)*x1.*x1;plot(x1,y1,'*')holdonplot(x1,y2,'r')

結(jié)果：

平方誤差：程序代碼：

symsabx=[2,3,4,7,8,10,11,14,16,18,19];fi=a+b./x;y1=[106.42,108.2,109.5,110,109.93,110.49,110.59,110.6,110.76,111,111.2];y=log(y1);fy=fi-y;fy2=fy.^2;J=sum(fy.^2);Ja=diff(J,a);Jb=diff(J,b);Ja1=simple(Ja);Jb1=simple(Jb);[a,b]=solve(Ja1,Jb1);f=exp(a).*exp(b./x);plot(x,y1,'*')holdonplot(x,f,'r')xlabel('x')ylabel('y')legend('數(shù)值點(diǎn)','擬合后的圖',0)>>s=0;fori=1:1:11s=s+(f(i)-y1(i)).^2;endvpa(s)結(jié)果平方誤差

擬合誤差平方和為：0.47193209,通過計(jì)算比較發(fā)現(xiàn)，通過對(duì)數(shù)變換之后得到的擬合結(jié)果比進(jìn)行二次函數(shù)擬合的結(jié)果要好很多。所以對(duì)于一些特殊的函數(shù)擬合，可以通過一些變換，進(jìn)行擬合。3.5練習(xí)思考1、整體插值有何局限性？如何避免？2、簡(jiǎn)述數(shù)據(jù)擬合與插值的異同。解：1.線性插值的局限性：整體插值的過程中，若有無(wú)效數(shù)據(jù)則整體插值后插值曲線的平方誤差會(huì)比較大，即在該數(shù)據(jù)附近插值曲線的震動(dòng)幅度較大。在插值處理前，應(yīng)對(duì)原始數(shù)據(jù)進(jìn)行一定的篩選，剔除無(wú)效數(shù)據(jù)。對(duì)于函數(shù)高次插值穩(wěn)定性差，而低次插值對(duì)于較大區(qū)間逼近精度又不夠，而且，隨著插值節(jié)點(diǎn)數(shù)目的增多，區(qū)域中間的確越來(lái)越逼近原函數(shù)，然而區(qū)域兩端發(fā)生了較大震蕩，即發(fā)生了Runge現(xiàn)象，這是由叉子余項(xiàng)的表達(dá)式所決定的。解決這一矛盾的有效方法就是采用分段低次代數(shù)插值。2.①不同點(diǎn)：用插值