數值分析第二次程序題

1、數值分析第二次程序題插值法1.對Runge函數R(x) -2在區(qū)間-1,1作以下插值逼近，并和R(x)的圖像進行比擬,125x并對結果進行分析。(1)以 x -1 ih,h 0.1,0 i 20為節(jié)點，Newton 插值圖1-0.7,0.7上的Newton 插值圖2-1,1上的Newton 插值由上圖可以看出，在區(qū)間-0.7Q.7上，插值多項式可以比擬好地逼近被插值函數。而當區(qū)間改為-1,1時，邊界附近插值多項式與被插值函數的差異很大。即出現了 Runge現象。由于邊界接近60的誤差，圖像中間局部的變化幾乎不可見。主要原因是被插值函數的任意階導數不能到達一致有界。其插值余項(n 1)R(x)(

2、n 1)!n 1(x)不趨近零。插值多項式不能收斂到被插值函數。牛頓差值函數x=-1:0.1:1;endfunction f=niudun(z,N,n)n=length(x);for t=1:nf=N(1,1);for i=1:nc(t)=N(t,t);x=-1:0.1:1;y(i)=1/(1+25*x(i)*x(i);endfor k=2:nendz=-1:0.001:1;a=1;N=zeros(n,n);m=size(z,2);for r=1:(k-1)N(:,1)=y'for i=1:ma=a*(z-x(r);for j=2:nRunge(i)=1/(1+25*z(i)*z(i)

3、;endfor k=j:nf(i)=niudun(z(i),N,n);f=f+N(k,k)*a;endendN(k,j)=(N(k,j-1)-N(k-1,j-1)/(x(k)-x( k-j+1);plot(z,Runge,'k' ,z,f,'r')主程序end如下列圖，使用 Chebyshev多項式零點構造的 Lagrange插值多項式比擬接近原函數，沒 有出現Runge現象，圖4為第一小問中的等距節(jié)點插值，可以明顯的看出以 Chebyshev多項式零點為插值點的優(yōu)勢。主要原因是其多項式誤差為間內一致收斂。f (X)-Ln(X)12n(n(n 1)，在區(qū)Lagr

4、ange 函數function lag=lagrange(z,x,y)for i=1:21l(i)=1;for j=1:21if j=il(i)=l(i)*(z -x(j)/(x(i) -x(j);endendend1=1：lag=y*l;主程序for i=1:21 x(22-i)=cos(2*i-1)*pi/42);endfor i=1:21 y(i)=1/(1+25*x(i)*x(i);endz=-1:0.001:1;m=length(z);for i=1:m f(i)=1/(1+25*z(i)*z(i); lag(i)=lagrange(z(i),x,y);endplot(z,f,

5、9;k',z,lag,'r')沒有出現Runge現象。但是可以 因為這一局部的函數的曲率較大，以xi-1 ih,h 0.1,0 i 20為節(jié)點，分段線性插值如以下列圖所示，分段線性插值多項式比擬接近原函數， 明顯地看到在區(qū)間-0.101中，線性插值的擬合度較低,也就是二階導數較大。由誤差估計公式F( x)h2max8局部的誤差較大圖5線性插值由上圖可以看出，三次樣條插值函數的曲線及其光滑，圖中并沒有將插值函數連起來， 否那么根本無法分辨出原函數和插值函數的圖像，說明得到的函數十分接近被插值函數。另外,題目要求自然樣條插值，也就是再兩端的二階導數為0，需在變成過程中加以注

6、意。x二1:0.1:1;n=len gth(x);for i=1: n y(i)=1/(1+25*x(i)*x(i);endfor i=1: n-1for i=1: n-2 u(i)=h(i)/(h(i+1)+h(i); r(i)=1 -u(i);endG=zeros(n-2, n-2);for i=1: n-2h(i)=x(i+1)-x(i); endG(i,i)=2;endfor i=2: n-2G(i,i-1)=u(i-1);G(i,i+1)=r(i-1);endd=zeros(1, n-2);for i=1: n-2d(i)=6*(y(i+2) -y(i+1)/h(i+1) -(y(i

7、+1)-y(i)/h(i)/ (h(i+1)+h(i);endd=d'M=Gd;M=0;M;0;for i=1: n-1z=x(i):0.01:x(i+1);m=le ngth(z);for j=1:ms(j)=M(i)*(x(i+1) -z(j)A3/0.6+M(i+1)*(z(j) -x(i)A 3/0.6+(y(i) -M(i)*0.01/6)*(x(i+1) -z(j)/0.1+(y(i+ 1)-M(i+1)*0.01/6)*(z(j) -x(i)/0.1;endplot(z,s,'* r'，MarkerSize',3)hold onend hold o

8、nz二1:0.01:1;for i=1:201f(i)=1/(1+25*z(i)*z(i); endplot(z,f, 'b')2.對函數:1 sillTTX-1 < <0,幾工)=COS7TX0 < x < 1/2:01/2 < x < 1在區(qū)間-1,1作以下插值逼近，并和被插值函數的圖像進行比擬，并對結果進行分析。(1)以 xi-1 ih,h 0.1,0 i 20 為節(jié)點，Newt on 插值首先對函數進行簡要分析，函數 f(x)是分段函數，并且在 x=0處不連續(xù)，對于插值計 算，只需要函數值，所以除了函數作圖和計算函數值有所不同以外，程

9、序的主體局部沒有明顯改動，所以將此題程序統(tǒng)一放在最后。本小題中圖8-,上的Newton 插值由上圖可以看出， 在區(qū)間卜0.7,0.7上，插值多項式可以已經無法較好地逼近被插值函數了，而當區(qū)間改為-1,1時，邊界附近插值多項式與被插值函數的差異迅速擴大。即出現了 Runge現象。由于邊界接近 1000的誤差，圖像中間局部的變化幾乎不可見。相比于第一 題Runge現象更為明顯。丄n1(x)可能無窮大。(n 1)!插值多項式不能收斂到被插值函數。Lagra nge 插值圖9 以Chebyshev多項式零點為插值點主要原因是被插值函數不連續(xù)，導致其插值余項 Rn(x)如下列圖，使用 Chebyshev

10、多項式零點構造的 Lagrange插值多項式比擬接近原函數，沒 有出現Runge現象，并且可以看出，在不連續(xù)點位置插值效果一般， 但是在函數兩端的擬合 效果明顯要好，說明使用 Chebyshev多項式零點構造的 Lagrange插值多項式在連續(xù)函數上 的應用效果更佳。以Xi-1 ih,h 0.1,0 i 20為節(jié)點，分段線性插值圖10 21個插值點線性插值圖11 201個插值點線性插值如以下列圖所示，分段線性插值多項式比擬接近原函數，沒有出現Runge現象。此例中我們看到了線性插值的強大優(yōu)勢，當原函數較為光滑，曲率較小，即使是分段函數對線性插值的影響也極為有限，當插值點個數擴大10倍到達201

11、個時，可以明顯的看出線性插值的優(yōu)勢所在。以為-1 ih,h 0.1,0 i 20為節(jié)點，三次自然樣條插值圖12三次自然樣條插值函數圖像由上圖可以看出，三次樣條插值函數的曲線及其光滑， 但是與其他多項式擬合一樣在不 連續(xù)點處存在較大的誤差，但是與第一二小問中的 Lagra nge插值多項式相比，三次樣條插 值可以更快的脫離不連續(xù)點的影響，并在其他位置上表現出很好的擬合效果。綜合以上2題我們可以初步得出這樣的結論：當函數連續(xù)光滑，使用Chebyshev多項式零點構造的 Lagrange插值多項式可以有效地 防止Runge現象，但三次樣條插值函數的曲線更為優(yōu)秀。但是當函數出現不連續(xù)點時，分 段線性插

12、值的優(yōu)勢明顯，可以在不連續(xù)段處到達很好的擬合效果，并且可以迅速脫離不連 續(xù)點的影響，所以在做函數插值時在斜率很大的局部可以考慮使用分段線性插值，其他部 分采用三次樣條效果最好。牛頓插值for j=2:nv=linspace(-1,0,100);x二1:0.1:1;for k=j:nu=sin(pi*v);n=length(x);plot(v,u,'k')for i=1:10N(k,j)=(N(k,j-1)-N(k-1,j-1)/(x(k)-x(hold ony(i)=sin(pi*x(i);k-j+1);v=linspace(0,0.5,50);endendu=cos(pi*v

13、);for i=11:15endplot(v,u,'k')y(i)=cos(pi*x(i);for t=1:nhold onendc(t)=N(t,t);v=linspace(0.5,1,50);for i=15:nendu=0;y(i)=0;z二0.1:0.01:0.1;plot(v,u,'k')endm=length(z);hold onfor i=1:mN=zeros(n,n);nd(i)=niudun(z(i),N,n);plot(z,nd, 'r')N(:,1)=y'end以Chebyshev多項式零點為插y(i)=cos(pi

14、*x(i);plot(v,u,'k')值點endhold onfor i=1:21endv=linspace(0,0.5,50);x(22-i)=cos(2*i-1)*pi/42);u=cos(pi*v);endz二1:0.001:1;plot(v,u,'k')for i=1:21m=length(z);hold onif x(i)<0for i=1:mv=linspace(0.5,1,50);y(i)=sin(pi*x(i);lag(i)=lagrange(z(i),x,y);u=0;elseifendplot(v,u,'k')y(i)=

15、0;v=linspace(-1,0,100);hold onelseu=sin(pi*v);plot(z,lag,'r')線性插值elseifendx二1:0.01:1;y(i)=0;z=-1:0.001:1;for i=1:201elsen=length(z);if x(i)<0y(i)=cos(pi*x(i);m=floor(z+1)/0.01)+1;y(i)=sin(pi*x(i);endfor i=1:n-1l(i)=y(m(i)+(y(m(i)+1) -y(m(i)/(x(u=sin(pi*v);hold onm(i)+1)-x(m(i)*(z(i) -x(m(

16、i);plot(v,u,'k')v=linspace(0.5,1,50);endhold onu=0;l(2001)=y(201);v=linspace(0,0.5,50);plot(v,u,'k')f(2001)=y(201);u=cos(pi*v);hold onv=linspace(-1,0,100);plot(v,u,'k')plot(z,l,'r')三次樣條插值for i=1:n-2s(j)=M(i)*(x(i+1) -z(j)A3/0.6+M(i+x=-1:0.1:1;G(i,i)=2;1)*(z(j) -x(i)A3

17、/0.6+(y(i) -M(i)*0.01n=length(x);end/6)*(x(i+1) -z(j)/0.1+(y(i+1)-M(i+1for i=1:21for i=2:n-2)*0.01/6)*(z(j) -x(i)/0.1;if x(i)<0G(i,i-1)=u(i-1);endy(i)=sin(pi*x(i);G(i,i+1)=r(i-1);plot(z,s,'* r'，MarkerSize',3)elseifendhold ony(i)=0;d=zeros(1,n-2);endelsefor i=1:n-2v=linspace(-1,0,100);y(i)=cos(pi*x(i);u=sin(pi*v);endd(i)=6*(y(i+2) -y(i+1)/h(i+1)-(y(i+plot(v,u,'k')end1)-y(i)/h(i)/(h(i+1)+h(i);hold onendv=linspace(