數(shù)值分析實(shí)驗(yàn)報(bào)告

1、 實(shí)驗(yàn)2.1 多項(xiàng)式插值的振蕩現(xiàn)象實(shí)驗(yàn)?zāi)康模涸谝粋€(gè)固定的區(qū)間上用插值逼近一個(gè)函數(shù)，顯然Lagrange插值中使用的節(jié)點(diǎn)越多，插值多項(xiàng)式的次數(shù)就越高。我們自然關(guān)心插值多項(xiàng)式的次數(shù)增加時(shí)，Ln(x)是否也更加靠近被逼近的函數(shù)。Runge給出的一個(gè)例子是極著名并富有啟發(fā)性的。實(shí)驗(yàn)容：設(shè)區(qū)間-1，1上函數(shù) f(x)=1/(1+25x2)。考慮區(qū)間-1,1的一個(gè)等距劃分，分點(diǎn)為 xi= -1 + 2i/n，i=0，1，2，n，則拉格朗日插值多項(xiàng)式為.其中，li(x)，i=0，1，2，n是n次Lagrange插值基函數(shù)。實(shí)驗(yàn)步驟與結(jié)果分析：實(shí)驗(yàn)源程序function Chap2Interpolation%

2、 數(shù)值實(shí)驗(yàn)二：“實(shí)驗(yàn)2.1：多項(xiàng)式插值的震蕩現(xiàn)象”% 輸入：函數(shù)式選擇，插值結(jié)點(diǎn)數(shù)% 輸出：擬合函數(shù)及原函數(shù)的圖形promps = '請(qǐng)選擇實(shí)驗(yàn)函數(shù)，若選f(x),請(qǐng)輸入f,若選h(x),請(qǐng)輸入h,若選g(x),請(qǐng)輸入g:'titles = 'charpt_2'result = inputdlg(promps,'charpt 2',1,'f');Nb_f = char(result);if(Nb_f = 'f' & Nb_f = 'h' & Nb_f = 'g')er

3、rordlg('實(shí)驗(yàn)函數(shù)選擇錯(cuò)誤！');return;endresult = inputdlg('請(qǐng)輸入插值結(jié)點(diǎn)數(shù)N:','charpt_2',1,'10');Nd = str2num(char(result);if(Nd <1)errordlg('結(jié)點(diǎn)輸入錯(cuò)誤！');return;endswitch Nb_f case 'f' f=inline('1./(1+25*x.2)'); a = -1;b = 1; case 'h' f=inline('x./

4、(1+x.4)'); a = -5; b = 5; case 'g' f=inline('atan(x)'); a = -5; b= 5;end x0 = linspace(a, b, Nd+1); y0 = feval(f, x0); x = a:0.1:b; y = Lagrange(x0, y0, x); fplot(f, a b, 'co'); hold on; plot(x, y, 'b-'); xlabel('x'); ylabel('y = f(x) o and y = Ln(x)-&

5、#39;);%-function y=Lagrange(x0, y0, x);n= length(x0); m=length(x);for i=1:m z=x(i); s=0.0; for k=1:n p=1.0; for j=1:n if(j = k) p = p*(z - x0(j)/(x0(k) - x0(j); end end s = s + p*y0(k); end y(i) = s;end實(shí)驗(yàn)結(jié)果分析(1)增大分點(diǎn)n=2，3，時(shí)，拉格朗日插值函數(shù)曲線如圖所示。 n=6 n=7 n=8 n=9 n=10從圖中可以看出，隨著n的增大，拉格朗日插值函數(shù)在x=0附近較好地逼近了原來的函數(shù)f

6、(x)，但是卻在兩端x= -1和x=1處出現(xiàn)了很大的振蕩現(xiàn)象。并且，仔細(xì)分析圖形，可以看出，當(dāng)n為奇數(shù)時(shí)，雖然有振蕩，但振蕩的幅度不算太大，n為偶數(shù)時(shí)，其振蕩幅度變得很大。通過思考分析，我認(rèn)為，可能的原因是f(x)本身是偶函數(shù)，如果n為奇數(shù)，那么Lagrange插值函數(shù)Ln(x)的最高次項(xiàng)xn-1是偶次冪，比較符合f(x)本身是偶函數(shù)的性質(zhì)；如果n為偶數(shù)，那么Lagrange插值函數(shù)Ln(x)的最高次項(xiàng)xn-1是奇次冪，與f(x)本身是偶函數(shù)的性質(zhì)相反，因此振蕩可能更劇烈。(2)將原來的f(x)換為其他函數(shù)如h(x)、g(x)，結(jié)果如圖所示。其中h(x), g(x)均定義在-5，5區(qū)間上，h(

7、x)=x/(1+x4)，g(x)=arctan x。h(x), n=7 h(x), n=8h(x), n=9 h(x), n=10g(x), n=7g(x), n=8g(x), n=9g(x), n=10分析兩個(gè)函數(shù)的插值圖形，可以看出：隨著n的增大，拉格朗日插值函數(shù)在x=0附近較好地逼近了原來的函數(shù)f(x)，但是卻在兩端x= -5和x=5處出現(xiàn)了很大的振蕩現(xiàn)象。并且，仔細(xì)分析圖形，可以看出，當(dāng)n為偶數(shù)時(shí)，雖然有振蕩，但振蕩的幅度不算太大，n為奇數(shù)時(shí)，其振蕩幅度變得很大。原因和上面f(x)的插值類似，h(x)、g(x)本身是奇函數(shù)，如果n為偶數(shù)，那么Lagrange插值函數(shù)Ln(x)的最高次項(xiàng)

8、xn-1是奇次冪，比較符合h(x)、g(x)本身是奇函數(shù)的性質(zhì)；如果n為奇數(shù)，那么Lagrange插值函數(shù)Ln(x)的最高次項(xiàng)xn-1是偶次冪，與h(x)、g(x)本身是奇函數(shù)的性質(zhì)相反，因此振蕩可能更劇烈。實(shí)驗(yàn)3.1 多項(xiàng)式最小二乘擬合實(shí)驗(yàn)?zāi)康模?編制以函數(shù)xkk=0,n；為基的多項(xiàng)式最小二乘擬合程序。實(shí)驗(yàn)容： 對(duì)表中的數(shù)據(jù)作三次多項(xiàng)式最小二乘擬合。xi-1.0-0.50.00.51.01.52.0yi-4.447-0.4520.5510.048-0.4470.5494.552取權(quán)函數(shù)wi1，求擬合曲線中的參數(shù)k、平方誤差2，并作離散據(jù)xi, yi的擬合函數(shù)的圖形。實(shí)驗(yàn)源程序function

9、 Chap3CurveFitting% 數(shù)值實(shí)驗(yàn)三:“實(shí)驗(yàn)3.1”% 輸出:原函數(shù)及求得的相應(yīng)插值多項(xiàng)式的函數(shù)的圖像以及參數(shù)alph和誤差rx0 = -1:0.5:2;y0 = -4.447 -0.452 0.551 0.048 -0.447 0.549 4.552;n = 3; % n為擬合階次alph = polyfit(x0, y0, n);y = polyval(alph, x0);r = (y0 -y)*(y0 -y)' % 平方誤差x = -1:0.01:2;y = polyval(alph, x);plot(x, y, 'k-');xlabel('

10、x'); ylabel('y0 * and polyfit.y-');hold onplot(x0, y0, '*')grid on;disp('平方誤差:', num2str(r)disp('參數(shù)alph:', num2str(alph)實(shí)驗(yàn)結(jié)果平方誤差:2.1762e-005參數(shù)alph:1.9991 -2.9977 -3.9683e-005 0.54912實(shí)驗(yàn)4.1 實(shí)驗(yàn)?zāi)康模簭?fù)化求積公式計(jì)算定積分. 實(shí)驗(yàn)題目：數(shù)值計(jì)算下列各式右端定積分的近似值. 實(shí)驗(yàn)要求： （1）若用復(fù)化梯形公式、復(fù)化Simpson公式和復(fù)化Ga

11、uss-Legendre I 型公式做計(jì)算，要求絕對(duì)誤差限為，分別利用它們的余項(xiàng)對(duì)每種算法做出步長(zhǎng)的事前估計(jì). （2）分別用復(fù)化梯形公式，復(fù)化Simpson 公式和復(fù)化Gauss-Legendre I 型公式作計(jì)算. （3）將計(jì)算結(jié)果與精確解做比較，并比較各種算法的計(jì)算量. 實(shí)驗(yàn)程序：1.事前估計(jì)的Matlab程序如下： （1）用復(fù)化梯形公式進(jìn)行事前估計(jì)的Matlab程序 format long g x=2:0.01:3; f=-4*(3*x.2+1)./(x.2-1).3; %二階導(dǎo)函數(shù) %plot(x,f) %畫出二階導(dǎo)函數(shù)圖像 x=2.0; %計(jì)算導(dǎo)函數(shù)最大值 f=-4*(3*x2+1)

12、/(x2-1)3; h2=0.5*10(-7)*12/f; h=sqrt(abs(h2) %步長(zhǎng) n=1/h; n=ceil(1/h)+1 %選取的點(diǎn)數(shù) format long g x=0:0.01:1; f=8.*(3*x.2-1)./(x.2+1).3; %二階導(dǎo)函數(shù) %plot(x,f) %畫出二階導(dǎo)函數(shù)圖像 x=1; %計(jì)算導(dǎo)函數(shù)最大值 f=8.*(3*x.2-1)./(x.2+1).3; h2=0.5*10(-7)*12/f; h=sqrt(abs(h2) %步長(zhǎng) n=1/h n=ceil(1/h)+1 %選取的點(diǎn)數(shù) format long g x=0:0.01:1; f=log(3

13、).*log(3).*3.x; %二階導(dǎo)函數(shù) %plot(x,f); %畫出二階導(dǎo)函數(shù)圖像 x=1; %計(jì)算導(dǎo)函數(shù)最大值 f=log(3)*log(3)*3x; h2=0.5*10(-7)*12/f; h=sqrt(abs(h2) %步長(zhǎng) n=1/h n=ceil(1/h)+1 %選取的點(diǎn)數(shù) format long g x=1:0.01:2; f=2.*exp(x)+x.*exp(x);%二階導(dǎo)函數(shù) %plot(x,f) %畫出二階導(dǎo)函數(shù)圖像 x=2; %計(jì)算導(dǎo)函數(shù)最大值 f=2.*exp(x)+x.*exp(x); h2=0.5*10(-7)*12/f; h=sqrt(abs(h2) %步長(zhǎng)

14、 n=1/h n=ceil(1/h)+1 %選取的點(diǎn)數(shù)估計(jì)結(jié)果步長(zhǎng)h及結(jié)點(diǎn)數(shù)n分別為 h = 0.8 n =1793 h =0.6 n =1827 h =0.9 n =2458 h =0.9 n =7020 （2）用復(fù)化simpson公式進(jìn)行事前估計(jì)的Matlab程序 format long g x=2:0.01:3; f=-2*(-72*x.2-24).*(x.2-1)-192*x.2.*(x.2+1)./(x.2-1).5;%四階導(dǎo)函數(shù) x=2.0; f=-2*(-72*x2-24)*(x2-1)-192*x2*(x2+1)/(x2-1)5; %計(jì)算導(dǎo)函數(shù)最大值 h4=0.5*10(-7)

15、*180*16/f; h=sqrt(sqrt(abs(h4) %步長(zhǎng) n=1/h; %求分段區(qū)間個(gè)數(shù) n=2*ceil(1/h)+1 %選取的點(diǎn)數(shù) format long g x=0:0.01:1; f=4*(-72*x.2+24).*(x.2+1)-192*x.2.*(-x.2+1)./(x.2+1).5;%四階導(dǎo)函數(shù) x=1; f=4*(-72*x2+24)*(x2+1)-192*x2*(-x2+1)/(x2+1)5; %計(jì)算導(dǎo)函數(shù)最大值 h4=0.5*10(-7)*180*16/f; h=sqrt(sqrt(abs(h4)%步長(zhǎng) n=1/h; %求分段區(qū)間個(gè)數(shù) n=2*ceil(1/h)

16、+1 %選取的點(diǎn)數(shù) format long g x=0:0.01:1; f=log(3)4*3.x;%四階導(dǎo)函數(shù) x=1; f=log(3)4*3.x;%計(jì)算導(dǎo)函數(shù)最大值 h4=0.5*10(-7)*180*16/f; h=sqrt(sqrt(abs(h4)%步長(zhǎng) n=1/h; %求分段區(qū)間個(gè)數(shù) n=2*ceil(1/h)+1 %選取的點(diǎn)數(shù) format long g x=1:0.01:2; f=4*exp(x)+x.*exp(x);%四階導(dǎo)函數(shù) plot(x,f) %畫出原函數(shù) x=2; f=4*exp(x)+x.*exp(x); %計(jì)算導(dǎo)函數(shù)最大值 h4=0.5*10(-7)*180*16

17、/f; h=sqrt(sqrt(abs(h4) n=1/h; %求分段區(qū)間個(gè)數(shù) n=2*ceil(1/h)+1 %選取的點(diǎn)數(shù) 估計(jì)結(jié)果步長(zhǎng)h及結(jié)點(diǎn)數(shù)n分別為 h =0.13411 n =47 h =0.76542 n =35 h =0.18433 n =29 h =0.18546 n =49 2. 積分計(jì)算的Matlab程序： format long g promps='請(qǐng)選擇積分公式，若用復(fù)化梯形，請(qǐng)輸入T，用復(fù)化simpson，輸入S，用復(fù)化Gauss_Legendre，輸入GL：' result=inputdlg(promps,'charpt 4',1,&

18、#39;T'); Nb=char(result); if(Nb='T'&Nb='S'&Nb='GL') errordlg('積分公式選擇錯(cuò)誤'); return; end result=inputdlg('請(qǐng)輸入積分式題號(hào)1-4：','實(shí)驗(yàn)4.1',1,'1'); Nb_f=str2num(char(result); if(Nb_f<1|Nb_f>4) errordlg('沒有該積分式'); return; end switch N

19、b_f case 1 fun=inline('-2./(x.2-1)');a=2;b=3; case 2 fun=inline('4./(x.2+1)');a=0;b=1; case 3 fun=inline('3.x');a=0;b=1; case 4 fun=inline('x.*exp(x)');a=1;b=2; end if(Nb='T')%用復(fù)化梯形公式 promps='請(qǐng)輸入用復(fù)化梯形公式應(yīng)取的步長(zhǎng)：' result=inputdlg(promps,'實(shí)驗(yàn)4.2',1,&

20、#39;0.01'); h=str2num(char(result); if(h<=0) errordlg('請(qǐng)輸入正確的步長(zhǎng)！'); return; end tic; N=floor(b-a)/h); detsum=0; for i=1:N-1 xk=a+i*h; detsum=detsum+fun(xk); end t=h*(fun(a)+fun(b)+2*detsum)/2; time=toc; t end if(Nb='S')%用復(fù)化Simpson公式 promps='請(qǐng)輸入用復(fù)化Simpson公式應(yīng)取的步長(zhǎng)：' resul

21、t=inputdlg(promps,'實(shí)驗(yàn)4.2',1,'0.01'); h=str2num(char(result); if(h<=0) errordlg('請(qǐng)輸入正確的步長(zhǎng)！'); return; end tic; N=floor(b-a)/h); detsum_1=0; detsum_2=0; for i=1:N-1 xk_1=a+i*h; detsum_1=detsum_1+fun(xk_1); end for i=1:N xk_2=a+h*(2*i-1)/2; detsum_2=detsum_2+fun(xk_2); end t=

22、h*(fun(a)+fun(b)+2*detsum_1+4*detsum_2)/6; time=toc; t end if(Nb='GL')%用復(fù)化Gauss_Legendre I %先根據(jù)復(fù)化Gauss_Legendre I公式的余項(xiàng)估計(jì)步長(zhǎng) promps='請(qǐng)輸入用復(fù)化Gauss_Legendre I 公式應(yīng)取的步長(zhǎng)：' result=inputdlg(promps,'實(shí)驗(yàn)4.2',1,'0.01'); h=str2num(char(result); if(h<=0) errordlg('請(qǐng)輸入正確的步長(zhǎng)!

23、9;); return; end tic; N=floor(b-a)/h);t=0; for k=0:N-1 xk=a+k*h+h/2; t=t+fun(xk-h/(2*sqrt(3)+fun(xk+h/(2*sqrt(3); end t=t*h/2; time=toc; t end switch Nb_f case 1 disp('精確解：ln2-ln3=-0.4054651081') disp('絕對(duì)誤差：',num2str(abs(t+0.4054651081); disp('運(yùn)行時(shí)間：',num2str(time); case 2 dis

24、p('精確解：pi=3.979') disp('絕對(duì)誤差：',num2str(abs(t-pi); disp('運(yùn)行時(shí)間：',num2str(time); case 3 disp('精確解：2/ln3=1.368') disp('絕對(duì)誤差：',num2str(abs(t-1.368); disp('運(yùn)行時(shí)間：',num2str(time); case 4 disp('精確解：e2=7.065') disp('絕對(duì)誤差：',num2str(abs(t-7.065); d

25、isp('運(yùn)行時(shí)間：',num2str(time); end1. 當(dāng)選用復(fù)化梯形公式時(shí)： （1）式運(yùn)行結(jié)果為： t =-0.351 精確解：ln2-ln3=-0.4054651081 絕對(duì)誤差：1.3944e-008 運(yùn)行時(shí)間：0.003 （2）式運(yùn)行結(jié)果為： t =3.336 精確解：pi=3.979 絕對(duì)誤差：3.9736e-008 運(yùn)行時(shí)間：0.005 （3）式運(yùn)行結(jié)果為： t = 1.861 精確解：2/ln3=1.368 絕對(duì)誤差：4.3655e-008 運(yùn)行時(shí)間：0.016 （4）式運(yùn)行結(jié)果為： t =7.610 精確解：e2=7.065 絕對(duì)誤差：2.0775e-

26、008 運(yùn)行時(shí)間：0.007 2. 當(dāng)選用復(fù)化Simpson公式進(jìn)行計(jì)算時(shí)： （1）式運(yùn)行結(jié)果為： t =-0.7519 精確解：ln2-ln3=-0.4054651081 絕對(duì)誤差：2.7519e-011 運(yùn)行時(shí)間：0.022 （2）式運(yùn)行結(jié)果為： t =3.979 精確解：pi=3.979 絕對(duì)誤差：0 運(yùn)行時(shí)間：0.021 （3）式運(yùn)行結(jié)果為： t =1.288 精確解：2/ln3=1.368 絕對(duì)誤差：9.2018e-012 運(yùn)行時(shí)間：0.019 （4）式運(yùn)行結(jié)果為： t =7.118 精確解：e2=7.065 絕對(duì)誤差：9.0528e-011 運(yùn)行時(shí)間：0.021 3. 當(dāng)選用復(fù)化G

27、auss-Legendre I型公式進(jìn)行計(jì)算時(shí)： （1）式運(yùn)行結(jié)果為： t =-0.5262 精確解：ln2-ln3=-0.4054651081 絕對(duì)誤差：4.7385e-012 運(yùn)行時(shí)間：0.023 （2）式運(yùn)行結(jié)果為： t =3.979 精確解：pi=3.979 絕對(duì)誤差：1.3323e-015 運(yùn)行時(shí)間：0.021 （3）式運(yùn)行結(jié)果為： t =1.754 精確解：2/ln3=1.368 絕對(duì)誤差：6.1431e-012 運(yùn)行時(shí)間：0.019 （4）式運(yùn)行結(jié)果為： t =7.046 精確解：e2=7.065 絕對(duì)誤差：1.441e-012 運(yùn)行時(shí)間：0.021 結(jié)果分析： 當(dāng)選用復(fù)化梯形公

28、式時(shí)，對(duì)步長(zhǎng)的事前估計(jì)所要求的步長(zhǎng)很小，選取的節(jié)點(diǎn)很多，誤差絕對(duì)限要達(dá)到時(shí)，對(duì)不同的函數(shù)n 的取值需達(dá)到1000-10000之間，計(jì)算量是很大。 用復(fù)化simpson公式對(duì)步長(zhǎng)的事前估計(jì)所要求的步長(zhǎng)相對(duì)大些，選取的節(jié)點(diǎn)較少，誤差絕對(duì)限要達(dá)到時(shí)，對(duì)不同的函數(shù)n的取值只需在10-100之間，計(jì)算量相對(duì)小了很多，可滿足用較少的節(jié)點(diǎn)達(dá)到較高的精度，比復(fù)化梯形公式的計(jì)算量小了很多。用復(fù)化simpson公式計(jì)算所得的結(jié)果比用復(fù)化梯形公式計(jì)算所得的結(jié)果精度高很多，而且計(jì)算量小。 當(dāng)選用Gauss-Lagrange I型公式進(jìn)行計(jì)算時(shí)，選用較少的節(jié)點(diǎn)就可以達(dá)到很高的精度。實(shí)驗(yàn)5.1 常微分方程性態(tài)和R-K法穩(wěn)

29、定性試驗(yàn)實(shí)驗(yàn)?zāi)康模嚎疾煜旅嫖⒎址匠逃叶隧?xiàng)中函數(shù)y前面的參數(shù)對(duì)方程性態(tài)的影響（它可使方程為好條件的或壞條件的）和研究計(jì)算步長(zhǎng)對(duì)R-K法計(jì)算穩(wěn)定性的影響。實(shí)驗(yàn)容及要求： 實(shí)驗(yàn)題目：常微分方程初值問題其中，。其精確解為。實(shí)驗(yàn)要求：本實(shí)驗(yàn)題都用4階經(jīng)典R-K法計(jì)算。（1）對(duì)參數(shù)a分別取4個(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）取參數(shù)a為一個(gè)絕對(duì)值不大的負(fù)值和兩個(gè)計(jì)算步長(zhǎng)，一個(gè)步長(zhǎng)使參數(shù)ah在經(jīng)典R-K法的穩(wěn)定域，另一個(gè)步長(zhǎng)在經(jīng)典R-K法的穩(wěn)定域外。

30、分別用經(jīng)典R-K法計(jì)算并比較計(jì)算結(jié)果。取全域等距的10個(gè)點(diǎn)上的計(jì)算值，列表說明。實(shí)驗(yàn)程序：Matlab程序如下：function charp5RK%數(shù)值試驗(yàn)5.1：常微分方程性態(tài)和R-K法穩(wěn)定性試驗(yàn)%輸入：參數(shù)a，步長(zhǎng)h%輸出：精確解和數(shù)值解圖形對(duì)比%clf;result=inputdlg('請(qǐng)輸入-50，50間的參數(shù)a:','實(shí)驗(yàn)5.1',1,'-40');a=str2num(char(result); if (a<-50|a>50) errordlg('請(qǐng)輸入正確的參數(shù)a!'); return;endresult=inputdlg('請(qǐng)輸入（0 1）之間的步長(zhǎng):','實(shí)驗(yàn)5.1'