插值法與數(shù)據(jù)擬合.doc

佛山科學技術(shù)學院實 驗 報 告課程名稱 數(shù)值分析 實驗項目 插值法 專業(yè)班級 姓名 學號 指導教師 成 績 日 期 5月12日 一、實驗目的1、學會Lagrange 插值、牛頓插值和 分段線性插值等基本插值方法；2、討論插值的Runge現(xiàn)象，掌握分段線性插值方法3、學會Matlab提供的插值函數(shù)的使用方法，會用這些函數(shù)解決實際問題。二、實驗原理1、拉格朗日插值多項式2、牛頓插值多項式3、分段線性插值三、實驗步驟1、用MATLAB編寫獨立的拉格朗日插值多項式函數(shù)；2、用MATLAB編寫獨立的牛頓插值多項式函數(shù)；3、利用編寫好的函數(shù)計算本章P66例1、P77例1的結(jié)果并比較；4、已知函數(shù)在下列各點的值為：0.20.40.60.81.00.980.920.810.640.38試用4次牛頓插值多項式對數(shù)據(jù)進行插值，根據(jù)，畫出圖形。5、在區(qū)間-1,1上分別取用兩組等距節(jié)點對龍格函數(shù)作多項式插值，對不同值，分別畫出插值函數(shù)及的圖形。6、下列數(shù)據(jù)點的插值01491625364964012345678(1)可以得到平方根函數(shù)的近似，在區(qū)間0,64上作圖。（2）用這9個點作8次多項式插值。四、實驗結(jié)果1、用MATLAB編寫獨立的拉格朗日插值多項式函數(shù)Lagrange 插值多項式源代碼function yi=Lagrange(x, y, xi)% Lagrange 插值多項式，其中% x - 向量，全部的插值節(jié)點% y - 向量，插值節(jié)點處的函數(shù)值% xi - 標量，自變量x% yi - xi 處的函數(shù)估計值n=length(x); m=length(y);if n=m error(The lengths of X and Y must be equal); return;endp=zeros(1,n);% 對向量p賦初值0for k=1:n t=ones(1,n); for j=1:n if j=k if abs(x(k)-x(j)eps error(the DATA is error!); return; end t(j)=(xi-x(j)/(x(k)-x(j); end end p(k)=prod(t);endyi=sum(y.*p);2、用MATLAB編寫獨立的牛頓插值多項式函數(shù)function yi=New_Int(x, y, xi)% Newton 基本插值公式，其中% x - 向量，全部的插值節(jié)點，按行輸入% y - 向量，插值節(jié)點處的函數(shù)值，按行輸入% xi - 標量，自變量x% yi - xi 處的函數(shù)估計值n=length(x); m=length(y);if n=m error(The lengths of X and Y must be equal); return;end% 計算均差表YY=zeros(n); Y(:,1)=y; % Y(:,1)表示矩陣中第一列的元素for k=1:n-1 for i=1:n-k if abs(x(i+k)-x(i) y=0.41075,0.57815,0.69675,0.88811; yi=New_int(x,y,0.596)yi =0.6319144055040004、已知函數(shù)在下列各點的值為：0.20.40.60.81.00.980.920.810.640.38試用4次牛頓插值多項式對數(shù)據(jù)進行插值，根據(jù)，畫出圖形。解： X=0.2:0.2:1.0; y=0.98,0.92,0.81,0.64,0.38;xx=0.2:0.08:1.0; m=length(xx); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,o);plot(xx,z,r*);hold off得到如下圖形：圖一 練習4的圖形5、在區(qū)間-1,1上分別取用兩組等距節(jié)點對龍格函數(shù)作多項式插值，對不同值，分別畫出插值函數(shù)及的圖形。解：a=-1;b=1;n=100;h=(b-a)/n; x=a:h:b;y=1./(1+25.*x.2); plot(x,y,k)其函數(shù)原圖形分別如下所示：圖二 龍格函數(shù)的圖形用龍格函數(shù)的Lagrange（）插值函數(shù)畫圖源程序當n =10時，有：function Runge(10)% Runge現(xiàn)象% n - 等距離節(jié)點a=-1; b=1; h=(b-a)/n;x=a:h:b; y=1./(1+25.*x.2);xx=a:0.01:b; yy=1./(1+25.*xx.2); m=length(xx); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,o);plot(xx,z,r-);hold off當n =20時，有：function Runge(10)% Runge現(xiàn)象% n - 等距離節(jié)點a=-1; b=1; h=(b-a)/n;x=a:h:b; y=1./(1+25.*x.2);xx=a:0.01:b; yy=1./(1+25.*xx.2); m=length(xx); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,o);plot(xx,z,r-);hold off其圖形分別如下所示： 圖三 Runge(10) 圖四 Runge(20)6、下列數(shù)據(jù)點的插值01491625364964012345678(1)可以得到平方根函數(shù)的近似，在區(qū)間0,64上作圖。（2） 用這9個點作8次多項式插值。function L8%平方根函數(shù)的8次多項式插值x=0 1 4 9 16 25 36 49 64; y=sqrt(x);m=length(x); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, x(i);endhold onplot(x,y,o); plot(x,y,o);xlabel(x);ylabel(y);plot(x,z,r-);hold off得到其圖形如下所示：圖五 的圖形x=0 1 4 9 16 25 36 49 64;y=0 1 2 3 4 5 6 7 8;xx=0:0.5:64; yy=sqrt(xx);m=length(xx);z=zeros(1,m);for i=1:mz(i)=Lagrange(x,y,xx(i);endplot(x,y,o,xx,yy,.,xx,z,r-);的8次Lagrange插值函數(shù)圖形圖六 的8次Lagrange插值函數(shù)圖形五、討論分析從下面這段代碼可以看出，利用mathlab自帶的增長型變量來作圖，是一種在c語言學習中沒有發(fā)現(xiàn)的方法。而且，代碼中還體現(xiàn)了一種比較作圖和數(shù)據(jù)共用的思維方法。例如：x=0 1 4 9 16 25 36 49 64;y=0 1 2 3 4 5 6 7 8;xx=0:0.5:64; yy=sqrt(xx);m=length(xx);z=zeros(1,m);for i=1:mz(i)=Lagrange(x,y,xx(i);endplot(x,y,o,xx,yy,.,xx,z,r-);利用以上簡單的程序，matble就可以做出圖形，方便觀察和研究。我們知道，利用計算機則可以方便地增加曲線上樣點的個數(shù)，隨著樣點的個數(shù)增加，讓學生觀察函數(shù)的圖像性質(zhì)，從而抽象出函數(shù)圖像的形狀。當然，作圖的話我們也可以利用zz超級畫板。下面我們完成這個過程：（1）在坐標系的屬性對話框中，選擇“畫坐標網(wǎng)格”和“顯示刻度”選項。（2）在作圖區(qū)空白處單擊鼠標右鍵，在快捷菜單中單擊【函數(shù)或參數(shù)方程曲線】命令，如下圖所示彈出函數(shù)作圖屬性對話框，在“”對應的編輯框中輸入：x，然后在下面的曲線屬性中，設置“曲線的點數(shù)”為：9，選擇曲線以“折線段”方式顯示，設置變量x的參數(shù)范圍為：0到64，選擇在折線段上“畫點”，設置點的大小為：2；最后單擊【確定】按鈕完成,同樣可以畫出函數(shù)圖像，通過改變曲線的點數(shù)，函數(shù)曲線越來光滑；六、改進實驗建議三次樣條插值直接用spline函數(shù)做。 邊界條件加在y的首尾，第一個表示y(x0)，最后一個表示y(xt)。 如果不加邊界條件，默認是not-a-knot邊界條件（注意不是自然邊界條件） 自然邊界條件的插值要用csape函數(shù)才能得到。 如果用interp1，則只能使用spline函數(shù)的默認邊界條件，即not-a-knot條件。 下面是例子 :x=0:3:9; y=x.*cos(x); xx=linspace(0,9); plot(x,y,o);%樣本點 hold on; plot(xx,interp1(x,y,xx,spline),r);%interp1只能使用默認邊界條件plot(xx,spline(x,0 y 0,xx),r:);%spline可以使用第一類邊界條件，這里y(0)=y(9)=0 pp=csape(x,y,second); plot(xx,fnval(pp,xx) %第二類邊界條件要用csape做，這里自然邊界條件 legend(樣本點,默認邊界條件,一階導數(shù)為0,自然邊界條件,location,south)hold off;另外，我們可以設： y=a0+a1*x+a2*x2+a(n-1)*x(n-1)，則其 n次拉格朗日函數(shù)插值程序：function LagrangesNs() %用于求過n點的拉格朗日n-1次插值多項式options=Name of data file;title=Lagranges_points;lineNo=2;def=Lagranges.dat;outval=inputdlg(options,title,lineNo,def);if isempty(outval)=1,return,endfilename=outval1;data=load(filename);x=data(:,1);y=data(:,2);lagrangesN(x,y);endfunction lagrangesN(x,y)%畫出已知n個點的位置plot(x,y,*);hold on%n次拉格朗日多項式為 y=a0+a1*x+a2*x2+a(n-1)*x(n-1) %其中a0 a1 a2a(n-1)為待求系數(shù)n=length(x);X=Vandermonde(x,1);A=Xy;%繪制插值函數(shù)圖象x1=linspace(0,max(x);x2=Vandermonde(x1,n);y1=x2*A;plot(x1,y1);hold on%顯示公式func=y= ,num2str(A(1);for i=2:n; b=+ ,num2str(A(i),*x,num2str(i-1); func=func,b;endtext(0.8,0.8,func);end%創(chuàng)建