MATLAB實(shí)現(xiàn)拉格朗日插值精編版

1、數(shù)值分析上機(jī)報(bào)告題目：插值法學(xué)號(hào)： 201014924姓名：靳會(huì)有1一、調(diào)用 MATLAB內(nèi)帶函數(shù)插值1、MATLAB內(nèi)帶插值函數(shù)列舉如下：interp1一維數(shù)據(jù)內(nèi)插 (查表法 )interpft使用 FFT 方法的一維數(shù)據(jù)內(nèi)插interp2二維數(shù)據(jù)內(nèi)插 (查表法 )interp3三維數(shù)據(jù)內(nèi)插 (查表法 )interpn多維數(shù)據(jù)內(nèi)插 (查表法 )spline三次樣條內(nèi)插meshgrid為三維繪圖產(chǎn)生X 和 Y 陣ndgrid為多維函數(shù)和內(nèi)插產(chǎn)生陣列g(shù)riddata數(shù)據(jù)網(wǎng)格2、取其中的一維數(shù)據(jù)內(nèi)插函數(shù)（interp1 ）為例，程序如下：其調(diào)用格式為：yi=interp1(x, y, xi)yi=

2、interp1(x, y, xi, method)舉例如下：x=0:10:100y=40 44 46 52 65 76 80 82 88 92 110;xi=0:1:100yi=interp1(x,y,xi,'spline')3、其他內(nèi)帶函數(shù)調(diào)用格式為：Interpft 函數(shù)：y=interpft(x,n)y=interpft(x,n,dim)interp2 函數(shù)：ZI=interp2(X, Y, Z, XI, YI)，ZI=imerp2(Z, ntimes)2ZI=interp2(Z, XI, YI)ZI=interp2(X, Y , Z, XI, YI, method)in

3、terp3VI=interp3(X,Y ,Z,V,XI,YI,ZI)VI=interp3(V , ntimes)VI=interp3(V ,XI,YI,ZI)VI=interp3(, method)InterpnVI=interpn(X1, X2, X3, Y1,V Y2, Y3,)VI=interpn(V , ntimes)VI=interpn(V , Yl, Y2, Y3,VI=interpn(), method)Splineyi=spline(x,y,xi)pp=spline(x,y)meshgridX,Y=meshgrid(x,y)X,Y=meshgrid(x)X,Y ,Z=meshg

4、rid(x,y,z)NdgridX1, X2, X3,=ndgrid(x1, x2, x3, )X1, X2, X3,=ndgrid(x)GriddataZI=griddata(x, y, z, XI, YI)XI, YI, ZI=griddata(x, y, z, xi, yi)=griddata(method)1Mfunction f = Language(x,y,x0)syms t l;if(length(x) = length(y)n = length(x);elsedisp('x 和 y 的維數(shù)不相等！ ');return;%檢錯(cuò)3endh=sym(0);for (i

5、=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ù)end在 MATLAB 中輸入：x=18 31 66 68 70 72 70;y=23 33 52 51 43 40 46;f=Language(x,y)pl

6、ot(x,y)結(jié)果為：f =Inf + (-t)*Inf - 54329.8*t2 + 1503.75*t3 - 22.2065*t4 + 0.16789*t5 - 0.000512106*t6圖形如下：4MATLAB 實(shí)現(xiàn)拉格朗日插值建立如下拉格朗日插值函數(shù)：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;end5畫圖程序如下：x=-

7、5:1:5;y=1./(1+x.2);x0=-5:0.001:5;y0=lagrange(x,y,x0);y1=1./(1+x0.2);plot(x0,y0,'r')hold onplot(x0,y1,'g')注：畫出的圖形為n =10 的圖形得到圖形如下：6n=10的圖像牛頓 K 次插值多項(xiàng)式一、實(shí)驗(yàn)?zāi)康模?、掌握牛頓插值法的基本思路和步驟。2、 培養(yǎng)編程與上機(jī)調(diào)試能力。二、牛頓插值法基本思路與計(jì)算步驟：給定插值點(diǎn)序列 （ xi , f ( xi ) , i0,1, n, 。構(gòu)造牛頓插值多項(xiàng)式N n (u) 。輸入要計(jì)算的函數(shù)點(diǎn) x, 并計(jì)算 N n (x)

8、的值，利用牛頓插值公式，當(dāng)增加一個(gè)節(jié)點(diǎn)時(shí)，只需在后面多計(jì)算一項(xiàng)，而前面的計(jì)算仍有用；另一方面N n ( x) 的各項(xiàng)系數(shù)恰好又是各階均差，而各階均差可用均差公式來計(jì)算。為的 一階均差。為的 k階均差。均差表：xk零 階 均一階均差二階均差三階均差7差Xf(X)00Xf(X)fX, X1101X2f(X 2)fX1,X 2fX0,X 1, X2Xf(X)fX, XfX, X,X fX ,X , X2332312301X3MMMMM牛頓插值法計(jì)算步驟：1 輸入 n 值及（ xi , f (xi ),i0,1, n, ；要計(jì)算的函數(shù)點(diǎn)x 。2 對(duì)給定的 x, 由N n ( x)f ( x0 ) (x

9、 x0 ) f x0 , x1( x x0 )( x x1) f x0 , x1, x2(x x0 )( x x1)( x xn 1 ) f x0 , x1, xn計(jì)算 N n ( x) 的值。3輸出 N n ( x) 。程序清單：functionc, d=newpoly(x, y)%牛頓插值的 MA TLAB 實(shí)現(xiàn)%這里 x 為 n 個(gè)節(jié)點(diǎn)的橫坐標(biāo)所組成的向量，y 為縱坐標(biāo)所組成的向量。%c 為所求的牛頓插值多項(xiàng)式的系數(shù)構(gòu)成的向量。n=length(x);% 取 x 的個(gè)數(shù)。d=zeros(n, n);% 構(gòu)造 nXn 的空數(shù)組。d(: , 1)=y'forj=2 : nfork=j

10、 : nd(k, j)=(d(k, j-1) - d(k-1, j-1) / (x(k)-x(k-j+1);endendc =d(n, n);for k=(n-1) : - 1 : 1c =conv(c, poly(x(k);% conv求積， poly(x) 將該多項(xiàng)式的系數(shù)賦給向量。m=length(c);c(m)=c(m)+d(k, k);end五、測(cè)試數(shù)據(jù)與結(jié)果：測(cè)試數(shù)據(jù)：（第三章習(xí)題第三題第2 題）f(x)=lnx的數(shù)值如表所示，構(gòu)造牛頓插值多項(xiàng)式并求ln0.53 的值。X0.40.50.60.70.8ln-0.916-0.693-0.51082-0.35776-0.223144x29114765解：由表可知 x0=0.4, x 1 =0.5, x 2=0.6, x 3=0.7, x 4=0.7，函數(shù)值：Y0=-0.916291,y1=-0.693147,y2=-0.510826, y3=-0.357765,y4=-0.2231448建立一個(gè)主程序np.mclcclearnewpoly(0.4,0.5,0.6,0.7,0.