東南大學(xué)數(shù)值分析上機(jī)作業(yè)

1、 數(shù)值分析上機(jī)報(bào)告院系：學(xué)號(hào)：姓名： 目錄作業(yè)1、舍入誤差與有效數(shù)11、函數(shù)文件cxdd.m12、函數(shù)文件cddx.m13、兩種方法有效位數(shù)對(duì)比14、心得2作業(yè)2、Newton迭代法21、通用程序函數(shù)文件22、局部收斂性3（1）最大值文件3（2）驗(yàn)證局部收斂性43、心得5作業(yè)3、列主元素Gauss消去法61、列主元Gauss消去法的通用程序62、解題中線性方程組73、心得8作業(yè)4、三次樣條插值函數(shù)81、第一型三次樣條插值函數(shù)通用程序：82、數(shù)據(jù)輸入及計(jì)算結(jié)果10作業(yè)1、舍入誤差與有效數(shù)設(shè)，其精確值為.（1） 編制按從小到大的順序，計(jì)算的通用程序；（2） 編制按從大到小的順序，計(jì)算的通用程序；（

2、3） 按兩種順序分別計(jì)算,并指出有效位數(shù)；（4） 通過本上機(jī)你明白了什么？程序：111、函數(shù)文件cxdd.mfunction S=cxdd(N)S=0;i=2.0;while(i<=N) S=S+1.0/(i*i-1); i=i+1;endscript運(yùn)行結(jié)果（省略>>）：S=cxdd(80)S= 0.7375772、函數(shù)文件cddx.mfunction S=cddx (N)S=0;for i=N:-1:2 S=S+1/(i*i-1);endscript運(yùn)行結(jié)果（省略>>）：S=cddx(80)S=0.7375773、兩種方法有效位數(shù)對(duì)比精確值函數(shù)：functio

3、n S=jqz(N)S=0.5*(1.5-1.0/N-1.0/(N+1);script運(yùn)行結(jié)果（省略>>）NS精確值從小到大從大到小值有效位數(shù)值有效位數(shù)1000.7400500.74005060.7400496100000.7499000.74990040.749852410000000.7499990.74999960.74985234、心得本題重點(diǎn)體現(xiàn)了數(shù)值計(jì)算中“大數(shù)吃小數(shù)”的問題，由于計(jì)算機(jī)計(jì)算的截?cái)嗵攸c(diǎn)，從大到小的計(jì)算會(huì)導(dǎo)致小數(shù)的有效數(shù)被忽略掉。從題中可以看出，看出按不同的順序計(jì)算的結(jié)果是不相同的，按從小到大的順序計(jì)算的值與精確值吻合，而按從大到小的順序計(jì)算的值與精確值有

4、較大的誤差。計(jì)算機(jī)在進(jìn)行數(shù)值計(jì)算時(shí)會(huì)出現(xiàn)“大數(shù)吃小數(shù)”的現(xiàn)象，導(dǎo)致計(jì)算結(jié)果的精度有所降低。作業(yè)2、Newton迭代法（1） 給定初值x0及容許誤差，編制Newton法解方程f(x)=0根的通用程序。（2） 給定方程f(x)=x3/3-x=0,易知其有三個(gè)根x1=，x2=0，x3=。由Newton方法的局部收斂性可知存在0，當(dāng)x0（，），Newton迭代序列收斂于根x2 ，試確定盡可能大的；試取若干個(gè)初始值，觀察當(dāng)x0（-，-1），（-1，），（，），（，1），（1，+）時(shí)，Newton序列是否收斂以及收斂于哪一個(gè)根。（3） 通過本上機(jī)題，你明白了什么？1、通用程序函數(shù)文件定義f(x)函數(shù)fun

5、ction f=fun(x)f=x3/3-x;end定義f(x)導(dǎo)函數(shù)function f=dfun(x)f=x*x-1;end定義求近似解函數(shù)function f,n=newton(x0,ep)flag=1;n=0;while(flag=1) x1=x0-fun(x0)/dfun(x0); n=n+1; if(abs(x1-x0)<=ep|n>100000) flag=0; end x0=x1;endf=x1;end script運(yùn)行結(jié)果clear; x0=input('請(qǐng)輸入初始值x0：'); ep=input('請(qǐng)輸入容許誤差：&#

6、39;); f,n=newton(x0,ep);fprintf('方程的一個(gè)近似解為:%fn',x1); 2、局部收斂性（1）最大值文件flag=1;k=1;x0=0;while flag=1 sigma=k*10-6; x0=sigma; k=k+1; m=0; flag1=1; while flag1=1 && m<=103 x1=x0-fun(x0)/dfun(x0); if abs(x1-x0)<10-6 flag1=0; end m=m+1; x0=x1; end if (flag1=1|abs(x0)>=10-6

7、) flag=0; endendfprintf('最大值為： %fn',sigma); 運(yùn)行結(jié)果為：最大值為：0.774597即得最大的為0.774597，Newton迭代序列收斂于根=0的最大區(qū)間為（-0.774597，0.774597）。（2）驗(yàn)證局部收斂性l 在x0（-，-1）區(qū)間，取以下初值，分別調(diào)用newton.m函數(shù)文件，得到結(jié)果如下：X0X1迭代次數(shù)-100-1.73205115-20-1.73205111-5-1.7320518-1.5-1.7320515結(jié)果顯示，以上初值迭代序列均收斂于-1.732051，即根。顯然，迭代格式初值的選擇對(duì)于迭代的收斂速度是至關(guān)

8、重要的，當(dāng)初值接近真實(shí)值的時(shí)候，迭代次數(shù)減少。l 在x0（-1，）區(qū)間，取以下初值，分別調(diào)用newton.m函數(shù)文件，得到結(jié)果如下：X0X1迭代次數(shù)-0.951.7320519-0.851.7320516-0.801.73205110-0.781.73205115計(jì)算結(jié)果顯示，迭代序列局部收斂于1.730251，即根。l 在x0（，）區(qū)間，取以下初值，分別調(diào)用newton.m函數(shù)文件，得到結(jié)果如下：X0X1迭代次數(shù)-0.700.0000005-0.200.0000003-0.050.00000030.050.00000030.200.00000030.700.0000005由newton1.m

9、的運(yùn)行過程表明，在整個(gè)區(qū)間上均收斂于0，即根。l 在x0（，1）區(qū)間，取以下初值，分別調(diào)用newton.m函數(shù)文件，得到結(jié)果如下：X0X1迭代次數(shù)0.80-1.732051100.90-1.73205170.95-1.73205190.98-1.73205112計(jì)算結(jié)果顯示，迭代序列局部收斂于-1.732051，即根。l 在x0（1，+）區(qū)間，取以下初值，分別調(diào)用newton.m函數(shù)文件，得到結(jié)果如下：X0X1迭代次數(shù)1.51.732051551.7320518201.732051111001.73205115結(jié)果顯示，以上初值迭代序列均收斂于1.732051，即根。綜上所述：(-，-1)區(qū)間

10、收斂于-1.73205，(-1，)區(qū)間局部收斂于1.73205，局部收斂于-1.73205，(-，)區(qū)間收斂于0，(，1)區(qū)間類似于(-1，)區(qū)間，(1,)收斂于1.73205。3、心得牛頓迭代法對(duì)于初值的選擇要求較高，因此，在牛頓迭代時(shí)可現(xiàn)通過簡單迭代法尋找相對(duì)準(zhǔn)確一些的值來進(jìn)行牛頓迭代。對(duì)于方程有多解的問題，Newton法求方程根時(shí)，牛頓迭代要考慮局部收斂的問題，迭代序列收斂于某一個(gè)根有一定的區(qū)間限制，在一個(gè)區(qū)間上，可能會(huì)局部收斂于不同的根。作業(yè)3、列主元素Gauss消去法對(duì)于某電路的分析，歸結(jié)為求解線性方程組RI=V。 32 -13 0 0 0 -10 0 0 0 -13 35 -9 0

11、 -11 0 0 0 0 0 -9 31 -10 0 0 0 0 0R= 0 0 0 -30 57 -7 0 -5 0 0 0 0 0 -7 47 -30 0 0 0 0 0 0 0 -30 41 0 0 0 0 0 0 -5 0 0 27 -2 0 0 0 -9 0 0 0 -2 29VT=-15，27，-23，0，-20,12，-7,7,10T（1） 編制解n階線性方程組Ax=b的列主元Gauss消去法的通用程序；（2） 用所編程序解線性方程組RI=V，并打印出解向量，保留5位有效數(shù)字；（3） 在本編程之中，你提高了那些編程能力。1、列主元Gauss消去法的通用程序函數(shù)：找每列的主元的函數(shù)

12、function B=zhuyuan(B,t,N,M)for i=0:N-1-t if B(N-i,t)>B(N-i-1,t) c=zeros(1,M); for j=1:M c(j)=B(N-i,j); B(N-i,j)=B(N-i-1,j); B(N-i-1,j)=c(j); end endend 進(jìn)行列消去的函數(shù)function B=xiaoqu(B,t,N,M)for i=t+1:N l=B(i,t)/B(t,t); for j=t:M B(i,j)=B(i,j)-l*B(t,j); endend進(jìn)行三角矩陣下的解函數(shù)function X=jie(X,B,N,M)for i=1:

13、N-1 s=B(N-i,M); for j=N-i+1:N s=s-B(N-i,j)*X(j); end X(N-i)=s/B(N-i,N-i);end執(zhí)行主程序：N=input('請(qǐng)輸入線性方程組的階數(shù): N='); M=input('請(qǐng)輸入增廣矩陣階數(shù): M='); b=zeros(1,N);A=zeros(N,N);A=input('請(qǐng)輸入系數(shù)矩陣:');b(1,:)=input('請(qǐng)輸入方程組的右端向量:'); b=b' B=A,b; for t=1:N-1 B=zhuyuan(B,t,N,M); B=xiaoqu

14、(B,t,N,M);end X=zeros(N,1); X(N)=B(N,M)/B(N,N); X=jie(X,B,N,M); X2、解題中線性方程組執(zhí)行程序，輸入矩陣A（即題中的矩陣R）和列向量b(即題中的V),如下：請(qǐng)輸入線性方程組的階數(shù): n=9請(qǐng)輸入增廣矩陣階數(shù): m=10請(qǐng)輸入系數(shù)矩陣A：A=31,-13,0,0,0,-10,0,0,0,-15;-13,35,-9,0,-11,0,0,0,0,27;0,-9,31,-10,0,0,0,0,0,-23;0,0,-10,79,-30,0,0,0,-9,0;0,0,0,-30,57,-7,0,-5,0,-20;0,0,0,0,-7,47,-

15、30,0,0,12;0,0,0,0,0,-30,41,0,0,-7;0,0,0,0,-5,0,0,27,-2,7;0,0,0,-9,0,0,0,-2,29,10;請(qǐng)輸入方程組的右端向量b：-15 27 -23 0 -20 12 -7 7 10得到如下結(jié)果：A = 31 -13 0 0 0 -10 0 0 0 -15 -13 35 -9 0 -11 0 0 0 0 27 0 -9 31 -10 0 0 0 0 0 -23 0 0 -10 79 -30 0 0 0 -9 0 0 0 0 -30 57 -7 0 -5 0 -20 0 0 0 0 -7 47 -30 0 0 12 0 0 0 0 0

16、-30 41 0 0 -7 0 0 0 0 -5 0 0 27 -2 7 0 0 0 -9 0 0 0 -2 29 10x = -0.28923 0.34544 -0.71281 -0.22061 -0.43040 0.15431 -0.057823 0.20105 0.290233、心得列主元Gauss最重要的就是如何通過找到最大主元，并交換行，如何進(jìn)行消去，這需要很細(xì)心的循環(huán)，和很復(fù)雜的嵌套，通過本題的編程，感覺個(gè)人對(duì)于這種多層嵌套循環(huán)的處理能力提高了很多。通過對(duì)軟件的編程，也更加理解了Gauss消去法的實(shí)質(zhì)。作業(yè)4、三次樣條插值函數(shù)（1） 編制求第一型3次樣條插值函數(shù)的通用程序；（2）

17、已知汽車門曲線型值點(diǎn)的數(shù)據(jù)如下：i012345678910xi012345678910yi2.513.34.044.75.225.545.785.45.575.75.8端點(diǎn)條件為=0.8，=0.2，用所編程序求車門的3次樣條差值函數(shù)S（x），并打印出S（i+0.5），i=0,1······9。1、第一型三次樣條插值函數(shù)通用程序：n=input('請(qǐng)輸入節(jié)點(diǎn)數(shù)n:');n=n+1;xn=zeros(1,n);yn=zeros(1,n);xn(1,:)=input('請(qǐng)輸入X的值：');yn(1,:)=in

18、put('請(qǐng)輸入Y的值：');dy0=input('請(qǐng)輸入邊界條件y(0):');dyn=input('請(qǐng)輸入邊界條件y(n):'); d=zeros(n,1);h=zeros(1,n-1);f1=zeros(1,n-1);f2=zeros(1,n-2);for i=1:n-1 h(i)=xn(i+1)-xn(i);%求一階差商 f1(i)=(yn(i+1)-yn(i)/h(i);endfor i=2:n-1 f2(i)=(f1(i)-f1(i-1)/(xn(i+1)-xn(i-1); %求二階差商 d(i)=6*f2(i);endd(1)=6*

19、(f1(1)-dy0)/h(1);d(n)=6*(dyn-f1(n-1)/h(n-1); A=zeros(n); %求M值u=zeros(1,n-2);r=zeros(1,n-2);for i=1:n-2 u(i)=h(i)/(h(i)+h(i+1); r(i)=1-u(i);endA(1,2)=1;A(n,n-1)=1;for i=1:n A(i,i)=2;endfor i=2:n-1 A(i,i-1)=u(i-1); A(i,i+1)=r(i-1);endM=Ad;syms xfor i=1:n-1 %求節(jié)點(diǎn)插值 Sx(i)=collect(yn(i)+(f1(i)-(M(i)/3+M(i

20、+1)/6)*h(i)*(x-xn(i)+M(i)/2*(x-xn(i)2+(M(i+1)-M(i)/(6*h(i)*(x-xn(i)3); Sx(i)=vpa(Sx(i),4);endS=zeros(1,n-1);for i=1:n-1 x=xn(i)+0.5; S(i)=yn(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i)*(x-xn(i) +M(i)/2*(x-xn(i)2+(M(i+1)-M(i)/(6*h(i)*(x-xn(i)3;enddisp('S(x)='); %結(jié)果輸出for i=1:n-1 fprintf(' %s (%d,%d)n

21、',char(Sx(i),xn(i),xn(i+1); disp('···········································

22、···');enddisp('S(i+0.5)');disp(' i x(i+0.5) S(i+0.5)')for i=1:n-1 fprintf(' %d %.4f %.4fn',i,xn(i)+0.5,S(i)end2、數(shù)據(jù)輸入及計(jì)算結(jié)果請(qǐng)輸入節(jié)點(diǎn)數(shù)n:10請(qǐng)輸入X的值：0 1 2 3 4 5 6 7 8 9 10請(qǐng)輸入Y的值：2.51 3.30 4.04 4.70 5.22 5.54 5.78 5.40 5.57 5.70 5.80請(qǐng)輸入邊界條件y(0):0.8請(qǐng)輸入邊界條件y(n):0.2S(x)=

23、- 0.008514*x3 - 0.001486*x2 + 0.8*x + 2.51 (0,1)··········································

24、··· - 0.0044579*x3 - 0.013654*x2 + 0.81217*x + 2.5059 (1,2)·····································

25、3;······· - 0.0036544*x3 - 0.018475*x2 + 0.82181*x + 2.4995 (2,3)·································&#

26、183;··········· - 0.040924*x3 + 0.31696*x2 - 0.18448*x + 3.5058 (3,4)·····························&#

27、183;··············· 0.10735*x3 - 1.4624*x2 + 6.9328*x - 5.9839 (4,5)··························&

28、#183;·················· - 0.26848*x3 + 4.1752*x2 - 21.255*x + 40.996 (5,6)······················

29、3;······················ 0.42659*x3 - 8.3361*x2 + 53.813*x - 109.14 (6,7)···················

30、83;························· - 0.26786*x3 + 6.2474*x2 - 48.272*x + 129.06 (7,8)················