數(shù)值分析實(shí)驗(yàn)題(華科)匯總

1、數(shù)值分析實(shí)驗(yàn)作業(yè)專業(yè)：姓名：學(xué)號(hào)：實(shí)驗(yàn) 2.1 多項(xiàng)式插值的振蕩現(xiàn)象問題提出 ：考慮在一個(gè)固定的區(qū)間上用插值逼近一個(gè)函數(shù)，顯然 Lagrange 插值中使用的節(jié)點(diǎn)越多，插值多項(xiàng)式的第 2 頁 共 13 頁次數(shù)就越高，我們自然關(guān)心插值多項(xiàng)的次數(shù)增加時(shí)，是極著名并富有啟發(fā)性的，設(shè)區(qū)間 -1,1 上函數(shù)f (x)實(shí)驗(yàn)內(nèi)容 ：考慮區(qū)間 -1 ，1的一個(gè)等距離劃分，分點(diǎn)為2ixi1n則拉格朗日插值多項(xiàng)式為L(zhǎng)n(x) 是否也更加靠近逼近的函數(shù)， Runge 給出的例子11 25x2i 0,1, 2,., nnLn(x)i011 25xi2li(x)其中， li (x) ， i=0,1,2, ,n 是 n

2、次 Lagrange 插值函數(shù)。實(shí)驗(yàn)要求 ：(1)選擇不斷增大的分點(diǎn)數(shù)目 n=2 ，3， 畫出原函數(shù) f(x)及插值多項(xiàng)式函數(shù) Ln(x) 在-1， 1上的圖像， 比較并分析實(shí)驗(yàn)結(jié)果。(2)選擇其他的函數(shù)，例如定義在區(qū)間-5， 5上的函數(shù)，xh(x) 4 ,g(x) arctan x1 x4重復(fù)上述的實(shí)驗(yàn)看其結(jié)果如何。解：以下的 f(x) 、h(x) 、 g(x) 的為插值點(diǎn)用“*”表示，朗格朗日擬合曲線用連續(xù)曲線表示。通過三個(gè)函數(shù)的拉格朗日擬合可以看到，隨著插值點(diǎn)的增加，產(chǎn)生(1) f(x)多項(xiàng)式求值的振蕩現(xiàn)象 n=2多項(xiàng)式求值的振蕩現(xiàn)象 n=4多項(xiàng)式求值的振蕩現(xiàn)象 n=6多項(xiàng)式求值的振蕩現(xiàn)

3、象 n=5多項(xiàng)式求值的振蕩現(xiàn)象 n=7y0.60.40.2-0.2-1f(x) lagrange(x)0.80.8 -0.6 -0.4 -0.21多項(xiàng)式求值的振蕩現(xiàn)象 n=8多項(xiàng)式求值的振蕩現(xiàn)象 n=9多項(xiàng)式求值的振蕩現(xiàn)象 n=10多項(xiàng)式求值的振蕩現(xiàn)象 n=11第 4 頁 共 13 頁多項(xiàng)式求值的振蕩現(xiàn)象 n=2多項(xiàng)式求值的振蕩現(xiàn)象 n=3(2) h(x)多項(xiàng)式求值的振蕩現(xiàn)象 n=4多項(xiàng)式求值的振蕩現(xiàn)象 n=5多項(xiàng)式求值的振蕩現(xiàn)象 n=7多項(xiàng)式求值的振蕩現(xiàn)象 n=6多項(xiàng)式求值的振蕩現(xiàn)象 n=8多項(xiàng)式求值的振蕩現(xiàn)象 n=9第 5 頁 共 13 頁多項(xiàng)式求值的振蕩現(xiàn)象 n=10多項(xiàng)式求值的振蕩現(xiàn)象

4、 n=111h(x)lagrange(x)432y01-1-2-3-4(3) g(x)多項(xiàng)式求值的振蕩現(xiàn)象 n=2多項(xiàng)式求值的振蕩現(xiàn)象 n=3多項(xiàng)式求值的振蕩現(xiàn)象 n=42多項(xiàng)式求值的振蕩現(xiàn)象 n=6多項(xiàng)式求值的振蕩現(xiàn)象 n=5y多項(xiàng)式求值的振蕩現(xiàn)象 n=7x第 7 頁 共 13 頁多項(xiàng)式求值的振蕩現(xiàn)象 n=821.50.5-0.5-1.5g(x) lagrange(x)-10.5g(x) lagrange(x)多項(xiàng)式求值的振蕩現(xiàn)象 n=9-0.5-1-1.5多項(xiàng)式求值的振蕩現(xiàn)象 n=10yx多項(xiàng)式求值的振蕩現(xiàn)象 n=11第 8 頁 共13 頁實(shí)驗(yàn) 3.1 最小二乘法擬合編制以函數(shù) xk kn

5、 0為基的多項(xiàng)式最小二乘擬合程序，并用于對(duì)表中的數(shù)據(jù)作三次多項(xiàng)式最小二乘擬合。xi-1.0-0.50.00.51.01.52.0yi-4.447-0.4520.5510.048-0.4470.5494.552n* * k 2取權(quán)數(shù) i 1,求擬合曲線 *k*xk 中的參數(shù) k ，平方誤差 2 ，并作離散數(shù)據(jù) xi,yik0的擬合函數(shù) y * (x) 的圖形。解：三次多項(xiàng)式的擬合曲線為：y (x) a0 a1x a2x2 a3x3此題中權(quán)函數(shù) (x) 1，即 W=(1,1,1,1,1,1,1)利用法方程 ATAa = ATY求解這個(gè)方程組，就可以得到系數(shù)a。解之得： 0 0.54912, 1 3

6、.9683 10 5, 2 2.9977, 3 1.9991 故擬合的函數(shù)為 :y 0.54912 3.9683 10 5x 2.9977x2 1.9991x3 , 平方誤差為： 2.176191667187105e-05擬合的函數(shù)圖像如下：-5離散值擬合曲線3次多項(xiàng)式擬合 ,平方誤差 =2.1762e-054321-1-2-3-4-0.5-10.5 xy0實(shí)驗(yàn) 5.1 常微分方程性態(tài)和 R-K 法穩(wěn)定性試驗(yàn)試驗(yàn)?zāi)康?： 考察下面的微分方程右端項(xiàng)中函數(shù) y 前面的參數(shù)對(duì)方程性態(tài)的影響 （它可使方程為好條 件的或壞條件的）和研究計(jì)算步長(zhǎng)對(duì) R-K 法計(jì)算穩(wěn)定性的影響。實(shí)驗(yàn)題目 ： 常微分方程初值

7、問題y y x 1, 0 x 1 y（0） 1 ,其中， 50 50 。其精確解為 y（x） e x x實(shí)驗(yàn)要求 ：（ 1）對(duì)于參數(shù) ，分別去四個(gè)不同的數(shù)值：一個(gè)大的正值，一個(gè)小的正值，一個(gè)絕對(duì) 值小的負(fù)值和一個(gè)絕對(duì)值大的負(fù)值。取步長(zhǎng) h 0.01 ，分別用經(jīng)典 R-K 法計(jì)算，將四組計(jì) 算結(jié)果畫在同一張圖上，進(jìn)行比較并說明相應(yīng)初值問題的性態(tài)。（2）對(duì)于參數(shù) 為一個(gè)絕對(duì)值不大的負(fù)值和兩個(gè)計(jì)算步，一個(gè)計(jì)算步使參數(shù)h 在經(jīng)典 R-K 法的穩(wěn)定域內(nèi)，另一個(gè)步長(zhǎng)在經(jīng)典的 R-K 法的穩(wěn)定域外。分別用經(jīng)典 R-K 法計(jì)算并 比較計(jì)算結(jié)果。取全域等距的 10 個(gè)點(diǎn)上的計(jì)算值，列表說明。解:對(duì)于 4階 R-

8、K法 絕對(duì)穩(wěn)定區(qū)為： 2.785 h 0這里 ，所以絕對(duì)穩(wěn)定區(qū)為： 2.785 h 0 （1）對(duì)于 h 0.01，絕對(duì)穩(wěn)定區(qū)：278.5 0a21-1-2h0.010.010.010.01微分方程數(shù)值解(2)對(duì)于20 ，穩(wěn)定區(qū) 0 h 0.1391a-20-20h0.010.158微分方程數(shù)值解 ,a=-20,h=0.011.1精確解 數(shù)值解微分方程數(shù)值解 ,a=-20,h=0.15765y432xy（精確解）數(shù)值解 y1 (a=-20,h=0.01 )y1-y數(shù)值解 y2( a=-20,h=0.15 )y1-y0.150.1997870.1997892.35E-061.5250001.325

9、2130.300.3024790.3024792.34E-072.1906251.8881460.450.4501230.4501231.75E-083.0496092.5994860.600.6000060.6000061.16E-094.1744633.5744570.750.7500000.7500007.23E-115.6648864.9148860.900.9000000.9000004.32E-127.6579696.757969可見 h=0.01 時(shí)，數(shù)值解穩(wěn)定 h=0.15 時(shí)，數(shù)值解不穩(wěn)定。第 15 頁 共 13 頁程序源代碼function testCharpt2_1%對(duì)數(shù)

10、值分析實(shí)驗(yàn)題第 2章第 1 題進(jìn)行分析promps=' 輸入 f 為選擇 f(x) ；輸入 h 為選擇 h(x) ；輸入 g 為選擇 g(x)' result=inputdlg(promps,' 請(qǐng)選擇實(shí)驗(yàn)函數(shù) ');chooseFunction=char(result);switch chooseFunctioncase 'f' f=inline('1./(1+25*x.2)');a=-1;b=1; nameFuc='f(x)'case 'h' f=inline('x./(1+x.4)&#

11、39;);a=-5;b=5 nameFuc='h(x)'case 'g' f=inline('atan(x)');a=-5;b=5 nameFuc='g(x)'end% promps2='n='% nNumble=inputdlg(promps2,' 請(qǐng)輸入分點(diǎn)數(shù) n'); nNumble=2:11for i=1:length(nNumble) x=linspace(a,b,nNumble(i)+1); y=feval(f,x);xx=a:0.1:b; yy=lagrange(x,y,xx) fig

12、urefplot(f,a,b,'*')hold onplot(xx,yy,'LineWidth',2)xlabel('x') ylabel('y') legend(nameFuc,'lagrange(x)') nameTitle=' 多項(xiàng)式求值的振蕩現(xiàn)象 ',' n=',num2str(nNumble(i) title(nameTitle,'FontSize',14);grid onendfunction yy=lagrange(x,y,xx)%s 實(shí)現(xiàn)拉格朗日插值%

13、輸入?yún)?shù) x ，y 分別為已知插值點(diǎn)的自變量和因變量%輸入?yún)?shù) xx 為擬合點(diǎn)的自變量值%輸出參數(shù) yy 為對(duì)應(yīng)自變量 xx 的擬合值 xLength=length(x);xxLength=length(xx);for i1=1:xxLengthyy(i1)=0;for i2=1:xLengthp=1;for i3=1:xLengthif(i2=i3) p=p*(xx(i1)-x(i3)/(x(i2)-x(i3);endendyy(i1)=yy(i1)+p*y(i2);endendfunction testCharpt3_1()%對(duì)數(shù)值分析實(shí)驗(yàn)題第 3章第 1 題進(jìn)行分析 %輸入?yún)?shù)：自變量

14、x，因變量 y %輸入?yún)?shù)：多項(xiàng)式擬合次數(shù) n clc clear format longx=-1.0,-0.5,0.0,0.5,1.0,1.5,2.0 y=-4.447,-0.452,0.551,0.048,-0.447,0.549,4.552 n=3A=;for i=1:length(x)A=A;1 x(i) x(i)2 x(i)3endA2=A'*A; a=inv(A2)*A'*y'% 多項(xiàng)式的系數(shù)% a=roundn(a,-6) yy=a(1)+a(2)*x+a(3)*x.2+a(4)*x.3; r=(y-yy)*(y-yy)' % 平方誤差clfhol

15、d on plot(x,y,'or'); x2=-1:0.01:2; y2=a(1)+a(2)*x2+a(3)*x2.2+a(4)*x2.3;plot(x2,y2,'LineWidth',2);legend('離散值 ','擬合曲線 ')xlabel('x');ylabel('y');title('3 次多項(xiàng)式擬合 ,平方誤差 =',num2str(r),'FontSize',14);grid onfunction testCharpt5_1%對(duì)數(shù)值分析實(shí)驗(yàn)題第 3章

16、第 1 題進(jìn)行分析 %輸入?yún)?shù)：參數(shù) a，步長(zhǎng) h %精確解和數(shù)值解圖形對(duì)比%第 1 問輸入a=2 1 -1 -2% 輸入 a 的取值h=0.01 0.01 0.01 0.01% 輸入 h 的取值%第 2 問輸入% a=-20 -20% 輸入 a 的取值% h=0.01 0.15% 輸入 h 的取值%func=inline('1+(y-x).*a');% 定義函數(shù)for i=1:length(a) x=0:h(i):1;% 求解區(qū)間 y=x;N=length(x); y(1)=1; for n=1:N-1k1=func(a(i),x(n),y(n); k2=func(a(i),x(n)+h(i)/2,y(n)+k1*h(i)/2); k3=func(a(i),x(n)+h(i)/2,y(n)+k2*h(i)/2); k4=func(a(i),x(n)+h(i),y(n)+k3*h(i) ; y(n+1)=y(n)+h(i)*(k1+2*k2+2*k3+k4)/6;% 數(shù)值解 endy0=exp(a