數(shù)值分析上級報告.doc

數(shù)值分析上機實驗報告函數(shù)插值方法課題五 函數(shù)插值方法一、 問題提出對于一元函數(shù)的n+1個節(jié)點值 j=(0,1,n)。試用Lagrange公式求其插值多項式或分段二次Lagrange插值多項式。數(shù)據(jù)如下：（1）0.40.550.650.800.951.050.410750.578150.696750.901.001.25382求五次Lagrange多項式，和分段三次插值多項式，計算,的值。（2）xj1234567yj0.3680.1350.0500.0180.0070.0020.001試構造Lagrange多項式，計算的值。0.165299 0.00213348二、要求 1、利用Lagrange插值公式編寫出插值多項式程序；2、給出插值多項式或分段三次插值多項式的表達式； 3、根據(jù)節(jié)點選取原則，對問題（2）用三點插值或二點插值，其結果如何； 4、對此插值問題用Newton插值多項式其結果如何。三、目的和意義 1、學會常用的插值方法，求函數(shù)的近似表達式，以解決其它實際問題； 2、明確插值多項式和分段插值多項式各自的優(yōu)缺點； 3、熟悉插值方法的程序編制； 4、如果繪出插值函數(shù)的曲線，觀察其光滑性。 四、實驗步驟1、根據(jù)（1）給出的實驗數(shù)據(jù)，編寫MATLAB程序，求出五次Lagrange多項式L5(x)，并計算,的值,畫出的函數(shù)圖。%五次lagrange 插值syms L5 Xx=0.4 0.55 0.65 0.80 0.95 1.05;y=0.41075 0.57815 0.69675 0.90 1.00 1.25382;L5=y(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)/(x(1)-x(2)*(x(1)-x(3)*(x(1)-x(4)*(x(1)-x(5)*(x(1)-x(6).+y(2)*(X-x(1)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)/(x(2)-x(1)*(x(2)-x(3)*(x(2)-x(4)*(x(2)-x(5)*(x(2)-x(6).+y(3)*(X-x(1)*(X-x(2)*(X-x(4)*(X-x(5)*(X-x(6)/(x(3)-x(1)*(x(3)-x(2)*(x(3)-x(4)*(x(3)-x(5)*(x(3)-x(6).+y(4)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(5)*(X-x(6)/(x(4)-x(1)*(x(4)-x(2)*(x(4)-x(3)*(x(4)-x(5)*(x(4)-x(6).+y(5)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(6)/(x(5)-x(1)*(x(5)-x(2)*(x(5)-x(3)*(x(5)-x(4)*(x(5)-x(6).+y(6)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)/(x(6)-x(1)*(x(6)-x(2)*(x(6)-x(3)*(x(6)-x(4)*(x(6)-x(5);L5x=vpa(L5)X=0.596;f1=subs(L5x)X=0.99;f2=subs(L5x)ezplot(L5x,0,2),title(五次lagrange插值函數(shù)),xlabel(X),ylabel(L5(X), hold on plot(x,y,r*) 計算結果為：L5x =385.79076923076923076923076923077*(X - 0.95)*(X - 0.4)*(X - 0.8)*(X - 0.65)*(X - 0.55) + 770.86666666666666666666666666667*(X - 0.95)*(X - 0.4)*(X - 0.8)*(X - 0.65)*(X - 1.05) - 1548.3333333333333333333333333333*(X - 0.95)*(X - 0.4)*(X - 0.8)*(X - 0.55)*(X - 1.05) + 1600.0*(X - 0.95)*(X - 0.4)*(X - 0.65)*(X - 0.55)*(X - 1.05) - 76.596736596736596736596736596737*(X - 0.95)*(X - 0.8)*(X - 0.65)*(X - 0.55)*(X - 1.05) - 1010.101010101010101010101010101*(X - 0.4)*(X - 0.8)*(X - 0.65)*(X - 0.55)*(X - 1.05)f1 =0.6257f2 =1.0542所以：f(0.596)= f1= 0.6257； f(0.99)= f2 =1.0542圖中紅色點為（1）中給出的實驗數(shù)據(jù)，從圖中可以看出，實驗數(shù)據(jù)點基本分布在插值函數(shù)的曲線上。2、對實驗數(shù)據(jù)（2），編寫MATLAB程序，求出六次Lagrange多項式L6(x)，并計算f(1.8),f(6.15)的值，畫出L6(x)的函數(shù)圖。%六次lagrange 插值syms L6 Xx=1 2 3 4 5 6 7;y=0.368 0.135 0.050 0.018 0.007 0.002 0.001;L6=y(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)*(X-x(7)/(x(1)-x(2)*(x(1)-x(3)*(x(1)-x(4)*(x(1)-x(5)*(x(1)-x(6)*(x(1)-x(7).+y(2)*(X-x(1)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)*(X-x(7)/(x(2)-x(1)*(x(2)-x(3)*(x(2)-x(4)*(x(2)-x(5)*(x(2)-x(6)*(x(2)-x(7).+y(3)*(X-x(1)*(X-x(2)*(X-x(4)*(X-x(5)*(X-x(6)*(X-x(7)/(x(3)-x(1)*(x(3)-x(2)*(x(3)-x(4)*(x(3)-x(5)*(x(3)-x(6)*(x(3)-x(7).+y(4)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(5)*(X-x(6)*(X-x(7)/(x(4)-x(1)*(x(4)-x(2)*(x(4)-x(3)*(x(4)-x(5)*(x(4)-x(6)*(x(4)-x(7).+y(5)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(6)*(X-x(7)/(x(5)-x(1)*(x(5)-x(2)*(x(5)-x(3)*(x(5)-x(4)*(x(5)-x(6)*(x(5)-x(7).+y(6)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(7)/(x(6)-x(1)*(x(6)-x(2)*(x(6)-x(3)*(x(6)-x(4)*(x(6)-x(5)*(x(6)-x(7).+y(7)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)/(x(7)-x(1)*(x(7)-x(2)*(x(7)-x(3)*(x(7)-x(4)*(x(7)-x(5)*(x(7)-x(6);L6x=vpa(L6)X=1.8;f1yp=subs(L6x)X=6.15;f2yp=subs(L6x)ezplot(L6x,0,7),title(六次lagrange插值函數(shù)),xlabel(X),ylabel(L6(X),hold onplot(x,y,r*)計算結果為：L6x = 0.00051111111111111111111111111111111*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 2.0)*(X - 4.0) - 0.001125*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 4.0) - 0.0005*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 2.0) + 0.00014583333333333333333333333333333*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 1.0)*(X - 2.0)*(X - 4.0) - 0.000016666666666666666666666666666667*(X - 7.0)*(X - 3.0)*(X - 5.0)*(X - 1.0)*(X - 2.0)*(X - 4.0) + 0.0010416666666666666666666666666667*(X - 7.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 2.0)*(X - 4.0) + 0.0000013888888888888888888888888888889*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 2.0)*(X - 4.0)f1yp = 0.1648f2yp =0.0013所以：f(1.8)= f1yp= 0.1648； f(6.15)= f2yp =0.0013圖中紅色點為（2）中給出的實驗數(shù)據(jù)分布圖。3、對于（2）給出的實驗數(shù)據(jù)，選取節(jié)點x2=2,y2=0.135； x4=4,y4=0.018；x6=6,y2=0.002；進行l(wèi)agrange插值，并畫出其函數(shù)圖，與六次lagrange進行比較。%三點Lagrange插值clearclcsyms L2 X1xyp=1 2 3 4 5 6 7;yyp=0.368 0.135 0.050 0.018 0.007 0.002 0.001;x=2 4 6;y=0.135 0.018 0.002;L2=y(1)*(X1-x(2)*(X1-x(3)/(x(1)-x(2)*(x(1)-x(3)+y(2)*(X1-x(1)*(X1-x(3)/(x(2)-x(1)*(x(2)-x(3).+y(3)*(X1-x(1)*(X1-x(2)/(x(3)-x(1)*(x(3)-x(2);L2x=vpa(L2)h=ezplot(L2,0,7);set(h,Color,k);set(h,linestyle,-);hold onplot(xyp,yyp,r*);hold on%六次Lagrange插值syms L6 Xx=1 2 3 4 5 6 7;y=0.368 0.135 0.050 0.018 0.007 0.002 0.001;L6=y(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)*(X-x(7)/(x(1)-x(2)*(x(1)-x(3)*(x(1)-x(4)*(x(1)-x(5)*(x(1)-x(6)*(x(1)-x(7).+y(2)*(X-x(1)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)*(X-x(7)/(x(2)-x(1)*(x(2)-x(3)*(x(2)-x(4)*(x(2)-x(5)*(x(2)-x(6)*(x(2)-x(7).+y(3)*(X-x(1)*(X-x(2)*(X-x(4)*(X-x(5)*(X-x(6)*(X-x(7)/(x(3)-x(1)*(x(3)-x(2)*(x(3)-x(4)*(x(3)-x(5)*(x(3)-x(6)*(x(3)-x(7).+y(4)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(5)*(X-x(6)*(X-x(7)/(x(4)-x(1)*(x(4)-x(2)*(x(4)-x(3)*(x(4)-x(5)*(x(4)-x(6)*(x(4)-x(7).+y(5)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(6)*(X-x(7)/(x(5)-x(1)*(x(5)-x(2)*(x(5)-x(3)*(x(5)-x(4)*(x(5)-x(6)*(x(5)-x(7).+y(6)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(7)/(x(6)-x(1)*(x(6)-x(2)*(x(6)-x(3)*(x(6)-x(4)*(x(6)-x(5)*(x(6)-x(7).+y(7)*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5)*(X-x(6)/(x(7)-x(1)*(x(7)-x(2)*(x(7)-x(3)*(x(7)-x(4)*(x(7)-x(5)*(x(7)-x(6);L6x=vpa(L6)ezplot(L6x,0,7),title(三點Lgrange插值和六次Lgrange插值的比較);%求函數(shù)值X1=1.8;X=1.8;k0=subs(L2x)k1=subs(L6x)X1=6.15;X=6.15;g0=subs(L2x)g1=subs(L6x)計算結果為：L2x =0.016875*(X1 - 6.0)*(X1 - 4.0) - 0.0045*(X1 - 6.0)*(X1 - 2.0) + 0.00025*(X1 - 2.0)*(X1 - 4.0)L6x =0.00051111111111111111111111111111111*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 2.0)*(X - 4.0) - 0.001125*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 4.0) - 0.0005*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 2.0) + 0.00014583333333333333333333333333333*(X - 7.0)*(X - 3.0)*(X - 6.0)*(X - 1.0)*(X - 2.0)*(X - 4.0) - 0.000016666666666666666666666666666667*(X - 7.0)*(X - 3.0)*(X - 5.0)*(X - 1.0)*(X - 2.0)*(X - 4.0) + 0.0010416666666666666666666666666667*(X - 7.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 2.0)*(X - 4.0) + 0.0000013888888888888888888888888888889*(X - 3.0)*(X - 6.0)*(X - 5.0)*(X - 1.0)*(X - 2.0)*(X - 4.0) X1 =1.8000X = 1.8000k0 = 0.1523k1 = 0.1648X1 = 6.1500X = 6.1500g0 = 0.0049g1 =0.0013所以：三點Lagrange插值：f(1.8)= k0= 0.1523； f(6.15)= k1=0.0049六次Lagrange插值：f(1.8)= g0= 0.1648； f(6.15)= g1=0.0013對比給出的數(shù)據(jù)：f(1.8)= 0.165299；f(6.15)= 0.00213348通過比較，六次lagrange能較好的擬合所給數(shù)據(jù)，取的數(shù)據(jù)點越多，插值函數(shù)越精確。4、對于（1）給出的實驗數(shù)據(jù)，編寫MATLAB程序，計算五次牛頓插值多項式p5(x)，畫出p5(x)的函數(shù)圖，計算p(0.596),p(0.99)值，并與lanrange插值法進行比較。%五次牛頓插值函數(shù)clearclcsyms p5 Xx=0.4 0.55 0.65 0.80 0.95 1.05;y=0.41075 0.57815 0.69675 0.90 1.00 1.25382;f0=(y(2)-y(1)/(x(2)-x(1);f01=(y(3)-y(2)/(x(3)-x(2);f02=(y(4)-y(3)/(x(4)-x(3);f03=(y(5)-y(4)/(x(5)-x(4);f04=(y(6)-y(5)/(x(6)-x(5);f1=(f01-f0)/(x(3)-x(1);f11=(f02-f01)/(x(4)-x(2);f12=(f03-f02)/(x(5)-x(3);f13=(f04-f03)/(x(6)-x(4);f2=(f11-f1)/(x(4)-x(1);f21=(f12-f11)/(x(5)-x(2);f22=(f13-f12)/(x(6)-x(3);f3=(f21-f2)/(x(5)-x(1);f31=(f22-f21)/(x(6)-x(2);f4=(f31-f3)/(x(6)-x(1);p5=y(1)+f0*(X-x(1)+f1*(X-x(1)*(X-x(2)+f2*(X-x(1)*(X-x(2)*(X-x(3)+f3*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)+f4*(X-x(1)*(X-x(2)*(X-x(3)*(X-x(4)*(X-x(5);p5x=vpa(p5)X=0.596,p1=subs(p5x),X=0.990,p2=subs(p5x),ezplot(p5x,0,2),title(),xlabel(X),ylabel(p5(X),hold on plot(x,y,r*)計算結果為：p5x =1.116*X + (X - 0.55)*(0.28*X - 0.112) + (0.99*X - 0.396)*(X - 0.65)*(X - 0.55) - 1.0*(15.30202020202020202020202020202*X - 6.1208080808080808080808080808081)*(X - 0.8)*(X - 0.65)*(X - 0.55) + (121.62635586635583706538454862311*X - 48.650542346542334826153819449246)*(X - 0.95)*(X - 0.8)*(X - 0.65)*(X - 0.55) - 0.03565X = 0.5960p1 = 0.6257X = 0.9900p2 = 1.0542所以：p6(0.596)=p1=0.6257， p6(0.99)= p2=1.05425、對于（2）給出的實驗數(shù)據(jù)，編寫MATLAB程序，計算六次牛頓插值多項式p6(x)，畫出p6(x)的函數(shù)圖，計算p(1.8),p(6.15)值，并與lanrange插值法進行比較。%六次牛頓插值syms p6 Xx=1 2 3 4 5 6 7;y=0.368 0.135 0.050 0.018 0.007 0.002 0.001;f0=(y(2)-y(1)/(x(2)-x(1);f01=(y(3)-y(2)/(x(3)-x(2);f02=(y(4)-y(3)/(x(4)-x(3);f03=(y(5)-y(4)/(x(5)-x(4);f04=(y(6)-y(5)/(x(6)-x(5);f05=(y(7)-y(6)/(x(7)-x(6);f1=(f01-f0)/(x(3)-x(1);f11=(f02-f01)/(x(4)-x(2);f12=(f03-f02)/(x(5)-x(3);f13=(f04-f03)/(x(6)-x(4);f14=(f05-f04)/(x(7)-x(5);f2=(f11-f1)/(x(4)-x(1);f21=(f12-f11)/(x(5)-x(2);f22=(f13-f12)/(x(6)-x(3);f23=(f14-f12)/(x(7)-x(4);f3=(f21-f2)/(x(5)-x(1);f31=(f22-f21)/(x(6)-x(2);f32=(f23-f22)/(x(7)-x(3);f4=(f31-f3)/(x(6)-x(1);f41=(f32-f31)/(x(7)-x(2);f5=(f41-f4)/(x(7)-x(1);p6=y(1)+f0*(X-x(1)+f1*(X-x(1)*(X-x(2)+